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ABSTRACT 


This publication contains the summaries for the Sixth Annual JPL Airborne Earth 
Science Workshop, held in Pasadena, California, on March 4r- 8 , 1996. The main 
workshop is divided into two smaller workshops as follows: 

•The Airborne Visible/Infrared Imaging Spectrometer (AVIRIS) workshop, on 
March 4-6. The summaries for this workshop appear in Volume 1. 

•The Airborne Synthetic Aperture Radar (AIRSAR) workshop, on March 6-8. 
The summaries for this workshop appear in Volume 2. 
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4X4 KENNADGB [K] MATRICES AND ITS INTERPRETATION IN TERMS OF THE 
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In polarimetry strict distinction must be made between optical polar ime try, 
defined for the forward propagation case in terms of the Jones/Mueller 
matrix calculus, versus radar polarimetry, defined for the (back) scattering 
case in terms of the Sinclair/Kennaugh matrix formulations. In this 
analysis, the optimal polarimetric contrast between a desirable (target) and 
undesirable (clutter) scattering feature is determined for radar polarimetry 
for which it is often desirable to select radar antenna polarization that 
maximize the contrast in received powers frcm such two distinct classes of 
scatterers. For example , a terrain-mapping POL-SAR imager may use 
polarization state manipulations to maximize the contrast between forested 
areas and farmland, whereas in a combined TOP-POL-SAR -*+■ POL-MTI imager, 
the polarization states may be chosen to optimize the "stealthy" aircraft 
target return relative to heterogeneous dense, dynamically changing 
background clutter. The fixed antenna polarization of such systems chosen 
for contrast optimization for either purpose may result in a target return 
that is small compared to external noise or internal receiver noise, clutter 
and speckle; and it may therefore be more appropriate to determine in a 
real-time image data take post-processing mode the polarization(s) that 
maximize(s) the signal- to-noise and its ratio with respect to speckle and 
clutter. A complete analysis of methods for the Sinclair [S] , the Graves 
[G] and the Kennaugh [K] matrices is provided for all of the scattering 
scenarios discussed here. Also, it is shwon why such a procedure applied to 

the covariance matrices [z 3 ' 4 ] for either the three(3)-dimensional 
reciprocal symmetric or the four (4)-dimensional bistatic (non- reciprocal) 
case in not useful; and that a group- theoretic expansion in terms of the 
"canonical eigenvalues" of the covariance matrix will instead provide the 
closely associated Polarimetric Entropy Coefficients (PEC H/A/a) of Cloude. 
Specifically, it is demonstrated that with the prior implementation of 
Cloude 's PEC (Polarimetric entropy H, polarimetric anisotropy A, derived 
frcm the eigenvalues of [Z] , and the polarimetric canonical target phase a 
derived frcm the associated eigenvectors) , it is possible to specify which 
class of OPCE coefficients is the roost effective one for truly optimal 
contrast enhancement. Various examples utilizing most recent image data 
take sets of the NAWC/ERIM-P3-UWB-T0PIF'E-P0LSAR and the NASA-JPL-AIRSAR , as 
well as the NASA/DARA/DASI-SIR-CA-SAR shuttle imaging systems will be 
presented for demonstrating the utility of the PEC and OPCE concepts in 
POL-SAR image anlaysis. 
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ON SOIL MOISTURE RETRIEVAL AND TARGET DECOMPOSITION 


Pascale C. Dubois(*,**) and Jakob van Zyl(*) 

(*): Jet Propulsion Laboratory, California Institute of Technology 
4800 Oak Grove Drive, Pasadena, CA 91109 

(**) Currently, visiting scientist with the Centre d'Etude des Environnements Terrestre et 

Planetaires, 78 Velizy, France 

E-mail: Pascale.dubois@cetp.ispl.fr Ph: 0 11 33 1 39 25 49 35 

1. INTRODUCTION 

In an earlier study, an empirical model was developed to infer soil moisture and surface 
roughness from radar data [1]. The inversion technique was extensively tested over bare surfaces 
by comparing the estimated soil moisture to in situ measurements. The overall RMS error in the 
soil moisture estimate was found to be 3.5 % and the RMS error in the RMS height estimate was 
less than 0.35 cm absolute for bare or slightly vegetated surfaces. However, inversion results 
indicate that significant amounts of vegetation cause the algorithm to underestimate soil moisture 
and overestimate RMS height. Among the areas over which the inversion cannot be applied, the 
areas with intermediate vegetation cover are of particular interest as both the vegetation and the 
underlying bare surface affect the backscatter. This paper concentrates mostly on these areas. 
Using the full polarimetric information and the Cloude target decomposition approach [2], three 
different components of the target backscattering can be isolated. One of these three components 
can be identified as the surface component in the case of intermediate vegetation cover. Once the 
surface component of the scattering is isolated, the bare surface inversion can then be applied. 

2. SOIL MOISTURE RETRIEVAL FOR BARE SURFACES 

The soil moisture retrieval method presented in [1] relies on the two following equations 
describing the hh-polarized and vv-polarized backscattering coefficients (7° hh and &’ vv as a 
function of 6 , the incidence angle, £ , the real part of the dielectric constant, h , the Root Mean 
square (RMS) height of the surface, k , the wave number and X , the wavelength in cm: 


M hh = lO- 2J5 ^^lO om£m9 (khsme) lA /l OJ 
sin 0 

j° w = IQ- 2 35 cos_gl 10° 0(^sin 0) 1 1 A 0 7 

sin 6 


( 1 ) 


The RMS height of the surface and the dielectric constant can easily be inverted from these 
two relations. Once the dielectric constant is known, the volumetric soil moisture can be 
computed using the Hallikainen curves [3] or the Brisco curves[4]. The inversion accuracy was 
extensively tested over a variety of sensors and sites as described in Figure 1. The overall RMS 
errors were found to be less than 0.4 cm in RMS height and 3.5% in soil moisture. The 
algorithm is optimized for bare surfaces and requires two co-polarized channels at a frequency 
between 1.5 GHz and 1 1 GHz. It gives best results for kh < 2.5 , jl v < 35% and 8 > 30° . 
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Figure 1: Radar derived soil moisture versus in situ measurements. 

Omitting the usually weaker hv-polarized returns makes the algorithm less sensitive to 
system cross-talk and system noise, simplifies the calibration process and adds robustness to the 
algorithm in the presence of vegetation. However, inversion results indicate that significant 
amounts of vegetation (NDVI > 0.4) cause the algorithm to underestimate soil moisture and 

overestimate RMS height. A simple criterion based on the <J° hv / cr^ ratio was developed to select 
the areas where the inversion is not impaired by the vegetation. In the following paragraphs, we 
will present a way to estimate the soil moisture for surfaces with intermediate vegetation cover. 


3. CLOUDE DECOMPOSITION 


Cloude showed that a general covariance matrix [T] can be decomposed as follows: 

[T] = Xik\ m ic] + ?i2k2*k2 + h3k3*kl + 7i4k4 9 kl (2) 

In (2), Xi , i= 1,2,3, 4 are the eigenvalues of the covariance matrix, , i= 1,2, 3,4 are the 

eigenvectors of [T], and k] means the adjoint (complex conjugate transposed) of In the 
monostatic case, for reciprocal media, the covariance matrix has one zero eigenvalue and the 
decomposition results in at most three covariance matrices on the right-hand side of (2). 

Also useful in our discussions later is Cloude's definition of target entropy. 
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( 3 ) 


4 2 

HT=Z 1 -Ptlog 4 (p,y where p, - 


As pointed out by Cloude, the target entropy is a measure of target disorder, with Ht = 1 for 
random targets and Hr = 0 for simple (single) targets. 


The covariance matrix for azimuthally symmetrical natural terrain in the monostatic case was 
shown [5] to follow the general form: 


M- 


yp 


o p 

Tj 0 

0 ? 


where 


J 


In n* \ {ShhSw) ^{ShvShv) , > ,_{SmS vv ) 

{ShhShh) ' p -W^)' 


(4) 

(5) 


The superscript * means complex conjugate, and all quantities are ensemble averages. The 
parameters c , £ , p and T] all depend on the size, shape and electrical properties of the 
scatterers, as well as their statistical angular distribution. It is easily shown that the eigenvalues 
of [T] are [6] 

A, = |(C + l + V(?-1) 2 + 4 PP*) (6) 

A 2 = ^(C + l-V(C-l) 2 + 4 pp* ) (7) 

A 3 = ci 7 (8) 

Because \T\ is a Hermitian matrix, the three eigenvalues are real and greater than zero. 


The corresponding three eigenvectors are 



(2p \ 

(C-i+Va) 2 

(Va+(C-d) 

0 

1 

v 7 

|(£-i + Va) 2 + 4 pp* 
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where A = (£-l) 2 + 4 pp . (12) 

The square root terms in front of the eigenvectors in (9-10) are there for normalization purposes 
and are always positive. We note that A is positive and that VA is always greater than |£ - 1| . 
Also note that we can write 

- = (13) 

foi (Va+(C-d) 

where D includes the square-root factors in front of the vectors in (9, 10) and is therefore always 
positive. It is easy to see then that the ratio in (13) is always real and negative. This means that 
the first two eigenvectors represent scattering matrices that can be interpreted in terms of odd and 
even numbers of reflections when Arg(p) is close to zero. This is the case for scattering 
dominated by the surface scattering term. 


4. SOIL MOISTURE RETRIEVAL FOR INTERMEDIATE VEGETATION COVER 


When surface scattering dominates, these conditions given by van Zyl [7] are met: 

Re(p) >$■ , !<1 and Arg(p) ~ 0 (14) 

The first eigenvalue, X\ corresponds then to the surface scattering term as fen is positive. It was 
pointed out in [1] that the presence of vegetation causes the inversion to underestimate the soil 
moisture and overestimate the roughness. In the following paragraph, we show that applying the 
inversion on the cross-sections corresponding to the first eigenvalue corrects this tendency by 
resulting in a higher soil moisture and a lower RMS height estimate than the straightforward 
inversion on the cross-sections corresponding to the original covariance matrix. 


From (14), we know that p is a positive real. We can then write (9) and (10) as: 


'01 s 

0 

_ _ 1 
J kl ~ r. 

(-P) 

0 

and fe 3 = 

V 

1 

^ 1 ) 

B 

v 1 , 


l°J 


where 

a = ■ t==~~ and /3 = ^ 

Va+(C-D Va-(C-i) 


(15) 


(16) 
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A, B, a and /? are real and positive. 


The original covariance matrix and the HH to VV ratio can be written as: 


[r] = 


A« 2 + A0 


2 q2 
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A 1 A 2 o 

+-TTP 


VA 




0 

A3 
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The covariance matrix corresponding to the first eigenvalue is: 


tn= 


%a 2 


0 

Ai 


-« 


0 

0 

0 


A, > 
— :a 
A 2 
0 

Ai 

A 2 


and - 


= a 


( 17 ) 


(18) 


(19) 


As long as £ is greater than 1 , a and the HH to VV ratio of the full covariance matrix is 
larger (closer to 1) than the HH to VV ratio of surface component of the covariance matrix. It 
follows that the soil moisture estimated from the first component of the decomposition will be 
higher than the soil moisture estimated from the full covariance matrix. It is also straightforward 

to see that the a° hh is larger for the original matrix than for the first component, resulting in a 
lower estimated value of the RMS height in the case of the surface component only. 

5. CONCLUDING REMARKS 


This paper describes a method to widen the domain of applicability of a soil moisture inversion 
algorithm previously published to include areas of intermediate vegetation cover. This method 
was tested on AIRSAR and SIR-C images and the results will be presented during the Workshop. 
In particular the following points will be clarified. 

< 7 ° 

The method applies to surfaces with an azimuthal symmetry, with < 1 and Arg(ShhSV) ~ 0 . 

<7w 

We will show that these conditions are met by most natural surfaces. 

In [1], the authors introduce the HV to VV ratio as a vegetation detector. We will evaluate how 
this criterion can be interpreted in the light of the Cloude decomposition and analyze the value of 
the decomposition step by identifying those surfaces rejected under the HV to VV ratio criterion 
for which the decomposition method allows an estimation of the soil moisture. 
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Abstract 

The NASA/JPL AIRSAR/TOPSAR instruments have the capability of collecting 
fully polarimetric radar data at three wavelengths (C, L and P-Bands) and dual antenna 
interferometry at C-Band, and more recently L-Band. In order to understand frequency 
and baseline dependent scattering effects in vegetated regions repeat pass interferometry 
data was collected for two vegetated regions in 1993. Portage Lake, Maine is a primarily 
coniferous forested region with some clear cutting from logging activities in the region. 
The second site at Innisfail, Australia borders a tropical rain forest and is situated 
adjacent to some major clear cut regions and banana plantations. Preliminary analysis of 
repeat pass data collected in these areas shows that the smaller the wavelength the greater 
the temporal decorrelation between passes, the longer the wavelength the greater the 
penetration depth for some types of vegetation canopy, yet for some vegetation canopy 
types, in particular for a banana plantation there appears to be no frequency dependent 
penetration into the canopy. 

Key words: SAR, interferometry, repeat-pass, vegetation 

Introduction 

Using two complex SAR images acquired from two spatially separated antennas 
collected simultaneously, single pass interferometry (SPI), or at different times, repeat 
pass interferometry (RPI), can be processed to yield elevation maps. The NASA/JPL (ri- 
band TOPSAR sensor is one example of an SPI interferometer [Zebker et al, 1992]. It is 
one component of the quadpolarization, three frequency airborne SAR (AIRSAR) 
operated by JPL. The processing of TOPSAR data into elevation maps is discussed in 
[Madsen et al, 1993]. 

It is well known that low frequency radars penetrate foliage better than radars at 
high frequencies, thus one might expect that interferometric measurements made at lower 
frequencies would yield height maps closer to the true ground surface. Until recently no 
multi-frequency SPI interferometric systems were available (TOPSAR has an L-Band 
interferometer added in 1995), hence this problem could only be addressed using RPI 
measurements. Moreover, to get interferometric measurements at wavelengths longer 
than L-band requires the use of RPI data. RPI can be used to generate height maps if 

1) the passes are acquired sufficiently close in time to avoid temporal decorrelation, 

2) the baseline vector, vector between the two antenna positions, is shorter than the 
critical baseline length, and 3) the baseline vector between the two passes is sufficiently 
constant. 


9 



Using repeat-pass interferometry data collected from the AIRSAR multi- 
frequency, multi-polarimetric radar provides a unique data set for studying frequency and 
polarization dependent phenomena on interferometric measurements in vegetated regions. 
Repeat pass data was collected near Portage Lake, Maine, a primarily coniferous forested 
area, and near Innisfail, Australia, a tropical forested area in Queensland. The radar 
operating parameters for the TOPSAR C-band mode, and the AIRSAR L and P-band 
polarimetric modes are given in Table I. 

Table I. AIRSAR/TOPSAR System Parameters 


Parameter 

C 

L 

P 

Wavelength (m) 

.0566 

.2422 

.6972 

Power (w) 

1000 

5000 

1000 

Bandwidth (MHz) 

40 

40 

40 

Pulse Spacing (m) 

.2572 

.5144 

.5144 

Antenna Length (m) 

1.5 

1.5 

1.5 

Baseline (m) 

2.51 

- 

- 

Altitude (m) 

8100 

8100 

8100 

Look Angle (°) 

20-60 

20-60 

20-60 

Range (km) 

8.5-17 

8.5-17 

8.5-17 

Ground Swath (km) 

12 

12 

12 


In order to process RPI data accurate platform position at each pulse event is 
crucial to determining the interferometric baseline as a function of time [Gray/Farris- 
Manning 1993], The NASA/JPL DC-8 is equipped with a Litton INU (subsequently 
upgraded to a Honeywell H764-G integrated GPS/INU), a TurboRogue differential GPS 
receiver, and pulse time tags derived from the GOES satellite. The INU data has an 
update rate of 50 Hz and gives accurate relative position measurements but tends to drift 
with time. The differential GPS is updated at a 1 Hz rate with estimated position accuracy 
of 10 cm. Precise platform position information is obtained by combining the INU data 
with the differential GPS data. These data are linked to the pulse data using the GOES 
time tag. 

Experiment Description 

Two sites were selected for study in 1993 using repeat-pass interferometry, one 
near Portage Lake in northern Maine and the other near Innisfail in Queensland, 
Australia. The first site has gently rolling terrain with elevations varying between 200 m 
and 300 m. The area is covered primarily with coniferous forest much of which is second 
growth due to widespread logging in the region. Vegetation heights vary from 1 m to 20 
m depending on the age and type of vegetation. There are many clear cut areas in varying 
stages of regrowth as well as open fields used for agriculture, mostly hay as the land is 
not particularly well suited for commercial farming. 

Considerable ground truth is available for this site since it has been used 
extensively by MIT Lincoln Laboratory for radar foliage penetration experiments in 
recent years. Measurements include tree height, canopy structure, biomass and dielectric 
for selected sites within a 20 km 2 area. In addition four 8 foot trihedral reflectors were 
deployed in open areas and surveyed using TurboRogue GPS receivers to be used as 
fiducial points in the images. 
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The NASA/JPL AIRSAR DC-8 flew eight passes over this site on June 24, 1993. 
Each pass was flown attempting to fly the same path each time in the hopes that at least 
one pair of passes would be suitable for RPI processing. As the onboard INU did not 
measure position and velocity accurately enough for RPI processing differential GPS 
measurements were acquired using a TurboRogue GPS receiver on the plane in 
conjunction with a reference receiver located near the site (and at the same location used 
as reference to survey the comer reflectors). Data was collected in the TOPSAR mode at 
C-band and in the AIRSAR 40 MHz quadpolarization mode at L and P-bands. 

The second site is located in Queensland, Australia and is situated along the 
boundary between a tropical rain forest and a clear cut region with cultivation of some 
crops such as bananas. The topography varies from relatively flat in open areas to very 
rugged in the tropical forest. The tropical canopy reaches heights of 30 m to 40 m with 
occasional clearings due to logging. 

Ground truth consists of 4 surveyed 8 foot trihedral reflectors deployed in open 
areas adjacent to the rain forest. In addition to the corner reflectors used as fiducial 
points, canopy measurements including tree heights and canopy structure were obtained 
by MIT Lincoln Laboratory. 

The NASA/JPL AIRSAR DC-8 flew eight passes over this site on September 4, 
1993. As described previously each pass was flown attempting to fly the same path each 
time in the hopes that at least one pair of passes would be suitable for RPI processing. 
Also as before TurboRogue differential GPS measurements were acquired using a 
TurboRogue GPS receiver on the plane in conjunction with a reference receiver located 
near the site (and at the same location used as reference to survey the corner reflectors). 
Data was collected in the TOPSAR mode at C-band and in the AIRSAR 40 MHz 
quadpolarization mode at L and P-bands. 

Data Processing 

Processing of the RPI data consists of four major steps, 1) the TurboRogue 
differential GPS measurements are combined with the INU data to obtain accurate 
platform position vectors, 2) the RPI passes over a given site are analyzed to determine 
which ones (if any) are suitable for interferometric processing, 3) the data are converted 
to an appropriate coordinate system and processed into interferograms, 4) residual motion 
measurements are made directly from the data and the scene is processed to a height map. 

TurboRogue data collected onboard the aircraft is processed at JPL using JPL 
derived ephemeris vectors for GPS satellites, and incorporates extensive models for 
ionospheric and tropospheric delays, yielding position vectors with an estimated accuracy 
of 10 cm [Webb et al, 1993], These data are then combined with the INU data by high 
pass filtering the INU data and combining it with GPS data using a least squares approach 
as described in [Madsen et al, 1993]. 

Next pairs of repeat pass lines are examined to find suitable pairs for RPI 
processing. One pass is arbitrarily chosen as reference and the platform position vectors 
are converted into a convenient coordinate frame, the "(s,c,h)" coordinate system where 
s, c and h are respectively the along track, cross track, and altitude coordinates of a "best 
fit" spherical approximation to the ellipsoid centered at the middle of the track. The 
(s,c,h) coordinate frame is not a Cartesian frame, but a spherical coordinate system 
designed to locally approximate the WGS-84 ellipsoid. This reference sphere is chosen to 
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be locally tangent to the ellipsoid at a reference point called the peg point, and to have 
radius equal to the radius of curvature in the along track direction as shown in Figure 1. 



Figure 1 . Approximating sphere used to define the (s,c,h) coordinate system. 

In the (s,c,h) system s denotes the distance along the reference curve from the peg point, 
which is a great circle in the approximating sphere, c is the distance along the arc 
perpendicular to the reference curve, and h is the height above the reference sphere. 

The transformation from radar mapping coordinates, (s,c,h), to WGS-84 
coordinates, (x,y,z) is described below. The radar mapping coordinates are defined 
relative to the sphere tangent to the ellipsoid at (0o,Ao), having radius r a which is the 
radius of curvature in the along track direction given by 


r e(^o) r n (A<)) 

r e (Ao) cos 2 (T|) + r n (Ao) sin 2 (rj) 


where q is the heading (actually track angle) and r e and r n are the radii of curvature in the 
east and north directions respectively given by 


r e (A) = - 

( 1 - e 2 sin 2 (a))? 


r n (^) = 


a (1 - e 2 ) 

( 1 - e 2 sin 2 (A.))? 


The point (00, Ao) will be referred to as the peg point and the sphere as the approximating 
sphere. Take the equator of the approximating sphere to be the reference curve (i.e. 
reference track), assumed to be a great circle on the approximating sphere. The prime 
meridian to the approximating sphere is the intersection of the approximating sphere and 
the plane determined by the unit normal vector, U, to the ellipsoid and the cross track 
vector C = UxT where T is the unit tangent vector to the reference curve. The radar 
mapping coordinates (s,c,h) are defined as the distance along the reference curve from the 
peg point, the distance from the reference curve along a meridian at a distance s from the 
peg point, and the height above the approximating sphere respectively. Let (x',y',z') be 
geocentric coordinates for the approximating sphere where the x' axis is U, the y' axis is 
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T and the z' axis is C. Then the transformation from mapping coordinates to x'y'z' 
coordinates is 


lx'} 


" (r a + h) cos(c^) cos(s 9 ) ’ 

y' 

= 

(r a + h) cos(c^) sin(s 0 ) 

\z'J 


> (r a + h) sin(c^) 


where c^ = c/r a and se = s/r a . The transformation from x'y'z' to xyz coordinates is an 
affine transformation of the form 


x 

y 

Z 


■ M xyz mFW . 

' iVi ENU M x y z 


X 

y 

\z / 


+ o 


where O is a translation vector, M£ is the transformation matrix from frame b to frame c, 
and the ENU .frame is a basis for the tangent space at the peg point with E a unit vector in 
the east direction and N a unit vector in the north direction. To obtain the transformation 
matrix from x'y'z' to ENU observe that they are related by a rotation about the U (or x') 
axis by the track angle r\. Thus the x'y'z' to ENU transformation matrix is 




0 

0 

1 


sin(ri) -cos(ri) 
cos(ri) sinCn) 
0 0 


The transformation matrix from ENU to xyz is just the matrix formed by the column 
vectors of E,N, and U in xyz coordinates. The transformation matrix is 


-sin(0o) -sin(Ao) cos(0o) cos(Xo) cos(0o) 

cos(0o) -sin(Ao) sin(0o) cos(Afr) sin(0o) 

0 cos(Ao) sin(Ao) 

The translation vector O is given by 

O = P-r a U 

where P=(r e (Xo)cos(A,o)cos(0o),r e (^O)cos(A,o)sin(0o),re(^O)(l-e 2 )sin(?to)) is the vector 
from the center of the ellipsoid to the peg point. 

After conversion of the platform position vectors to (s,c,h) coordinates of the 
reference path, the data are aligned in time and the instantaneous baseline vector can be 

computed. The instantaneous baseline vector is computed as follows. For each point on 

the reference track the instantaneous imaging plane is computed from the yaw and pitch 
angles. The intersection of this plane with the second track gives one candidate baseline 
vector. A second candidate results from taking the imaging plane from the point and 
finding its intersection on the first track. Finally, the intersection of the imaging plane at 
this point with the second track is computed. The reference track is given by the average 
of these four vectors and the normal to the imaging plane determined by the average of 
the normals to the four imaging planes. This is illustrated in Figure 2. The instantaneous 
baseline is the difference between the intersections of the imaging plane containing a 
point on the reference track with each of the two tracks in the RPI pair. Baselines for each 


M xyz 

m ENU 
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pair of passes are examined for length and convergence rates. RP1 pairs with baseline 
lengths from 15-50 meters and small convergence rates are the most desired. 


Track 1 Track 2 



Figure 2. Determining the reference path and instantaneous baseline. 

After selection of an appropriate RPI pair the data are motion compensated first to 
dual reference tracks and then to a single reference track as described in [Stevens et al, 
March 1995]. Using the single reference track approach the offsets between the two 
single look complex (SLC) images used to form the interferogram should be zero. Height 
reconstruction uses a modified form of the height reconstruction equations outlined in 
[Madsen et al, 1993]. 

Results 

At Portage Lake, Maine one pair of passes taken 30 minutes apart was suitable 
for RPI processing and the baseline as a function of along track position is shown in 
Figure 3. Amplitude images of all three passes are shown in Figure 4. Note that the flat 
open areas have progressively less backscatter as the wavelength increases whereas in the 
vegetated regions all three frequencies are relatively bright. Figure 5 shows C, L and P- 
band RPI inteferograms. One cycle of color corresponds to 2n radians of phase. The flat 
earth fringes are removed during the motion compensation process except for fringes due 
to topography and residual baseline errors. Figure 6 shows the TOPSAR generated 
topographic map. Note the close agreement of the P-band fringes with the TOPSAR 
height map. The higher azimuth fringing visible in both the C and L-band interferograms 
is a result of residual baseline errors. The magnitude of residual baseline can be estimated 
from the along track offsets between the two single look complex images used to form 
the interferogram. These offsets would be zero if knowledge of the baseline were perfect. 
A plot of the azimuth offsets as a function of along track position for five ranges is shown 
in Figure 7. Using a simple sinusoidal model for cross track position error as a function of 
along track position shows that the magnitude of the residual baseline error is 
approximately 6 cm. A subsequent paper will discuss estimation of these errors from the 
data. 


A comparison of the interferometric correlation for the three frequencies is quite 
interesting and is shown in Figure 8. In this figure blue denotes a correlation of .25 or 
below and yellow a correlation 1 .0. Note the almost complete decorrelation of the C-band 
interferogram in the vegetated areas while the clear cut areas are highly correlated. This 
can be compared with the TOPSAR correlation map shown in Figure 6 where the 
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correlation is generally above .9 in all areas. There is a dramatic increase in correlation in 
the vegetated regions for L and P-bands, however correlation in the flat open areas is 
quite poor particularly at P-band as a result of the low backscatter. We conclude that the 
dramatic drop in decorrelation as wavelength decreases is due to temporal decorrelation. 
Small motions in the leaves and branches on the order of one tenth of a wavelength cause 
the C-band to decorrelate faster than either L or P-bands. 

Because of the relatively short baseline of 15 m, the fringe frequency is quite 
small and unwrapping the P-band data is not too difficult, however because the baseline 
is almost completely horizontal the ambiguity height (height change for one cycle of 
phase) varies from 60 m in the near range to 600 m in the far range. This greatly reduces 
the height acuity in the far range. Figure 9 shows a P-band RPI height map. The height 
accuracy varies from 3-30 m from near to far range making height comparisons with the 
TOPS AR DEM useless in the far range. 

For the Innisfail, Australia data one pass was found to be suitable for RPI 
processing, however the baseline in this case is approximately 100 m. Figure 10 shows 
the TOPSAR generated DEM and L and P-band RPI interferograms. The TOPSAR DEM 
is shown in the upper left corner with one cycle of color equal to 150 meters, the 
ambiguity height at L-band. No attempt was made to remove residual baseline errors 
from the interferograms, hence there are slight ramps in both range and azimuth. The 
bright region in the upper left portion of the TOPSAR DEM is a banana plantation. 
Figure 1 la shows a cut through the trees, from which a tree height of approximately 7-8 
meters is inferred. After removing a tilt from the unwrapped P-band interferogram a cut 
through the banana plantation is shown in Figure 1 lb. Using an ambiguity height of 50 m 
for P-band the phase difference corresponds to roughly 8 m, in good agreement with the 
C-band data. Figure 10 shows the difference between the L and P band interferograms 
where the P-band phase has been scaled by the ratio of the wavelengths to match the L- 
band fringe rate. A cut through the phase difference plot is shown in Figure 1 lc. The 
phase difference corresponds to less than a meter. This indicates there is no differential 
penetration as a function of wavelength for the banana grove. 

Conclusions 

Preliminary analysis of repeat pass data collected in these areas shows that the 
smaller the wavelength the greater the temporal decorrelation between passes, the longer 
the wavelength the greater the penetration depth for some types of vegetation canopy, yet 
for some vegetation canopy types, in particular for a banana plantation there appears to be 
no frequency dependent penetration into the canopy for the three frequencies. 
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1. INTRODUCTION 

It has been demonstrated and recognized that radar interferometry is a promising method for 
the determination of digital elevation information (Zebker et al., 1994; Madsen et al., 1995) 
and terrain slope (Wegm tiller et al., 1993) from Synthetic Aperture Radar (SAR) data. An 
important application of Interferometric SAR (InSAR) data in areas with topographic varia- 
tions is that the derived elevation and slope can be directly used for the absolute radiometric 
calibration of the amplitude SAR data (van Zyl et al., 1993; Holecz et al., 1994), as well as 
for scattering mechanisms analysis (van Zyl, 1993; Holecz et al., 1995). On the other hand 
polarimetric SAR data has long been recognized as permitting a more complete inference of 
natural surfaces than a single channel radar system. In fact, imaging polarimetry provides the 
measurement of the amplitude and relative phase of all transmit and receive polarizations. 

On board the NASA DC-8 aircraft, NASA/JPL operates the multifrequency (P, L and C bands) 
multipolarimetric radar AIRSAR (van Zyl et al., 1992). The TOPSAR, a special mode of the 
AIRSAR system, is able to collect single-pass interferometric C- and/or L-band VV polarized 
data (van Zyl et al., 1995). A possible configuration of the AIRSAR/TOPSAR system is to 
acquire single-pass interferometric data at C-band VV polarization and polarimetric radar data 
at the two other lower frequencies. The advantage of this system configuration is to get digital 
topography information at the same time the radar data is collected. The digital elevation 
information can therefore be used to correctly calibrate the SAR data. This step is directly 
included in the new AIRSAR Integrated Processor. This processor uses a modification of the 
full motion compensation algorithm described by Madsen et al. (1993). However, the Digital 
Elevation Model (DEM) with the additional products such as local incidence angle map, and 
the SAR data are in a geometry which is not convenient, since especially DEMs must be referred 
to a specific cartographic reference system. Furthermore, geocoding of SAR data is important 
for multisensor and/or multitemporal purposes. In this paper, a procedure to geocode the new 
AIRSAR/TOPSAR data is presented. As an example an AIRSAR/TOPSAR image acquired 
in 1994 is geocoded and evaluated in terms of geometric accuracy. 
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2. THE AIRSAR INTEGRATED PROCESSOR 

A complete description of the AIRSAR/TOPSAR products is given in van Zyl (1995). Here we 
briefly present the different output data types. The Integrated AIRSAR Processor produces 
two data products: the standard AIRSAR products in frame and synoptic format, and the 
integrated TOPSAR product. The last one consists of a number of different data types, namely: 

• Digital Elevation Model , which represents the elevation of the terrain above a spher- 
ical approximation to the World Geodetic System 84 (WGS84). The data are stored as 
two bytes. 

• C-band VV polarized calibrated amplitude data as acquired with the top antenna in 
the TOPSAR mode. The radiometric corrections are performed taking into account the 
topography when removing the antenna patterns and the scattering areas. The data are 
stored as two bytes. 

• Local incidence angle map derived using the digital elevation model information. The 
data are stored as one byte. 

• Correlation map contains the normalized correlation coefficient between the two C- 
band interferometric channels. The data are stored as one byte. 

• L-band polarimetry data contain the calibrated AIRSAR compressed Stokes matrix 
data at L-band. Each pixel is represented by ten bytes. 

• P-band polarimetry data contain the calibrated AIRSAR compressed Stokes matrix 
data at P-band. Each pixel is represented by ten bytes. 

The selected radar mapping coordinate system (s, c, h r ) is defined relative to the sphere tangent 
to the WGS84 ellipsoid at longitude X p and latitude (p p - referred to as the peg point and 
approximately corresponding to the center of the image - having a radius r a , which is the 
radius of curvature in the along track direction (Hensley, 1993). In this system, s and c are 
distances in along and across direction, while h r is the height above the approximating sphere. 


3. FROM THE RADAR MAPPING TO MAP COORDINATE SYSTEM 

This section describes the transformation from the radar mapping coordinates ( s,c,h r ) into 
the map coordinates (x,y, h 0 ). An accurate transformation of the derived DEM data into map 
coordinates is only possible by using a forward transform. In fact, by applying a backward 
transform, which has the advantage of fewer number of calculations, the pixel location in the 
map coordinate system can be inaccurate, because the pixel altitude h 0 is usually unknown. 
Thus, the following geodetic and cartographic transformation steps must be carried out: 

1. From radar mapping ( s,c,h r ) into global Cartesian coordinates (X,Y,Z) 

2. From global Cartesian into local Cartesian coordinates (X', Y', Z') 

3. From local Cartesian into geographic coordinates (\,<p,h e ) 

4. From geographic into map coordinates ( x,y,h 0 ) 
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A forward transform implies that the output product can contain holes, if inappropriately 
implemented. To avoid gaps in the output data, four map coordinates, corresponding to four 
DEM corners in the radar mapping coordinates, define the orientation and the surface of the 
pixel. The area is then filled using a polygon fill algorithm. This procedure is subsequently 
applied to each pixel of the input data. In order to be computationally efficient, the transfor- 
mation of radar coordinate positions for each pixel to the map coordinate system is stored in 
a two dimensional look-up table. Thus, all additional input products, such as local incidence 
angle map, correlation map, and SAR data can be very easily and efficiently geocoded using 
the same look-up table. 

The first step is to transform the DEM data from the radar coordinates into the global Cartesian 
coordinates (in this case the WGS84 system), namely: 


where 
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r e and r n are the East and North radius of curvature at the peg point, 77 is the heading (or 
track angle), c v = c/r a , s\ = s/r a , O is a translation vector, Mi and M 2 are transformation 
matrices, a and e 2 are the semi-major axis and the ellipticity of the WGS84 ellipsoid (a = 
6378137 m, e 2 = 0.00669437999015). Basically, to specify the transformation from the radar 
mapping coordinates to the Earth centered Cartesian coordinate system, only three parameters 
are required, namely the latitude and longitude of the peg point, and the heading of the 
reference curve at the peg point. 
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Since the ellipsoid origin of a given cartographic reference system is not necessarily coincident 
with the Earth center and their shape and size may vary considerably, geodetic datums must 
be considered. The datum transformation is described by the three dimensional Helmert 
transformation (Frei et aL, 1993), namely: 


( x ) 

i Aa: ^ 


( X'\ 

Y 

= Ay 

+ (1 + m)D 

Y' 

\ z ) 

\A* J 


\Z' ) 


where 


cos / 3 cos 7 cos a sin 7 + sin a sin (3 sin 7 sin a sin 7 — cos a sin (3 cos 7 
D=| — cos (3 sin 7 cos a cos 7 — sin a sin (3 sin 7 sin a cos 7 + cos a sin f3 sin 7 
sin f3 — sin a cos ( 3 cos a cos / 3 


and (Ax,Ay,Az) are translation parameters, m a scaling parameter, and (a,/3, 7 ) three rota- 
tion parameters. 

Now, the transformation from the local Cartesian coordinates into the geographic coordinates 
can be carried out, namely: 


<p = arctan 
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Z' + (e') 2 b sin 3 9 
d — e 2 acos 3 0 
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a 2 -9 
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(1 — e 2 sin 2 (p) 2 

d = [(x') 2 + (f') 2 ]' 


9 = arctan 
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Note that a and b are the semi-major and semi-minor axes of the ellipsoid used in a given coun- 
try. Furthermore, not all countries provide the solution of the 7 parameters (i.e. Ax,Ay,Az, 
m,a,l3, 7 ). In those cases, usually, only the three translation parameters must be known. 
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The last step is to transform the geographic coordinates into a desired map coordinate system. 
The much used projections are the Universal Transfer Mercator, the Oblique Mercator, the 
Lambert Conformal Conic Projection, and the Stereographic Projection. A detailed description 
is given in Snyder (1987). As an example, we consider here the UTM system: 

x — x 0 T k\A A T k^AX^ -f- fc 5 AA 5 
y = yo T &o T &2AA 2 -j- &4AA 4 T &gAA 6 


where x a and y 0 are the false origin coordinates, AA = (A — A 0 ) with X 0 the geographic longitude 
of the central meridian, and ko...ke are coefficients (see Snyder, 1987). Note, x refers to the 
Easting coordinates, while y to the Northing coordinates. It should be noted that to obtain 
orthometric heights h a , ellipsoidal heights h e need to be corrected for the geoidal heights h g 
according to h 0 = h e — h g . 


4. RESULTS 


The SAR data used in this study were collected by the NASA/JPL AIRSAR/TOPSAR system 
in 40 MHz mode on April 12, 1994 over the area of Atlanta airport (33.76°N, 83.59°W), Georgia. 
The data have a postprocessing pixel spacing of 10 meters and cover an area of around 100 
km 2 . The DEM was geocoded to the UTM coordinate system zone 16 (central meridian 87°) 
taking into account the North American Datum 1927 (Clarke ellipsoid 1866, Ax = -9 m with 
cr Ax = 5 m, Ay = 161 m with <7 Ay = 5 m, Az = 179 m with oaz — 8 m). 


Figure 1 represents the location accuracy for 8 selected ground control points without us- 
ing (left) and by using some independent ground control points (right). Figure 2 shows the 
geocoded DEM (a) with the corresponding correlation map (b), the derived local incidence an- 
gle map (c), the C-band VV-polarized amplitude data (d), the L-band HH-polarized amplitude 
data (e), and the L-band HV-poIarized amplitude data (f ). 
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Figure 1: Positioning accuracy without (left) and with ground control points (right). 
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e) f) 

Figure 2: Geocoded AIRSAR/TOPSAR products in UTM zone 16 a) DEM, b) Correlation 
map, c) Local incidence angle map, d) C-VV data, e) L-HH data, f) L-HV data. 
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An important prerequisite for accurate geocoding is precise flight path data. Currently, the 
NASA DC-8 aircraft is equipped with the Honeywell H-764G Inertial Navigation System (INS) 
with an embedded Global Positioning System (GPS) receiver. This system, which has an 
acquisition rate of 50 Hz, meets a position accuracy of < 16 meters, a pitch and roll accuracy 
of < 0.01 degrees, and a heading accuracy of < 0.02 degrees (Honeywell, 1993). Note that this 
GPS performance is obtained during a period of normal GPS coverage with a Position Dilution 
of Precision of 3.2 or less, and P(Y)-Code receiver mode. In other words, using this flight 
path data accuracy, and also considering the uncertainties of the datums’ shift parameters, 
a geocoding precision of ±1 pixel is expected, assuming a postprocessing pixel spacing of 10 
meters. 

Figure 1 (left) clearly shows in the y— direction, i.e. N-S direction, an average location error 
of around 80 meters (corresponding to 8 pixels in the DEM data), while in the a:— direction 
it is about 1 pixel. However, by using 1 ground control point, these inaccuracies are reduced, 
resulting in a final positioning precision on the order of ±1 pixel (Figure 1 right). A possible 
explanation of the observed error in y — direction could be figured in a time synchronization 
problem between the GPS measurements and the SAR data. In fact, the flight track angle rj in 
this case was of 0°, i.e. S-N direction, and therefore positioning inaccuracies in the y — direction 
were mainly affected. It should be noted that the position of the 8 ground control points was 
determined by using a 1:24,000 scale map of the U.S. Geological Survey. Therefore, some minor 
inaccuracies due to the map scale must also be considered. 


5. CONCLUSIONS 

In this paper we presented a procedure to geocode the new AIRSAR/TOPSAR data. Fur- 
thermore, it was shown that geocoding of these data is well defined, without using, as is often 
done, transformations based on polynomials, which can strongly affect the geometric accuracy. 
Using geodetic and cartographic transforms, the geocoding can basically be carried out with- 
out the need of ground control points, which sometimes can be difficult, especially for those 
areas where the availability of topographic maps is almost impossible. At this time, we are 
testing this procedure for data sets acquired during the 1994 and 1995 campaigns in order to 
figure out the above mentioned error. In addition, we are planning to optimize the geocoding 
procedure without using the intermediary (s,c, h r ) radar mapping system. This will allow a 
straight geocoding from slant range to a given reference map coordinate system, in order to 
deliver to the users geocoded AIRSAR/TOPSAR data. 
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1. INTRODUCTION 

A new technique has been developed to measure ocean surface slopes by polarimetric SAR. This 
technique has the potential of estimating surface slope in the along-track direction at an accuracy of a 
fraction of a degree. Experimental results using NASA/JPL AIRSAR P-, L-, C-band data from the Gulf 
Stream’90 Experiment showed a sudden slope change at a converging current front. 

The capability of SAR for the imaging of ocean waves, currents, surfactants, slicks, etc. has been 
studied in great detail ( Hasselman et al, 1985, Alpers et al, 1981). However, few studies have been 
devoted to the understanding of the polarimetric aspect of the SAR response to waves, currents and other 
ocean surface features. The polarimetric SAR simultaneously collects the backscattered responses from 
HH, HV, VH and VV polarizations, from which a Mueller matrix can be formed. From the Mueller 
matrix, the inherent scattering mechanism can be revealed by the polarization signature (van Zyl, 1987). 
Schuler et al (1993) have shown that, based on the stochastic surface tilts suggested by Cox and Munk 
(1954), the polarimetric SAR response of an ambient ocean surface agrees with the tilted Bragg model ( 
Valenzuela, 1968). They also showed that the co-polarization signature of a current front displayed a 
higher degree of depolarization induced by the composite nature of scattering mechanisms, such as tilted- 
Bragg and specular, and the scatter caused by breaking waves. 

For idealized Bragg scattering, the peak value of the co-polarization signature plot occurs at a 
polarization orientation angle of 90°, and an ellipticity angle of 0°(i.e., at VV polarization). However, a 
tilted surface in the azimuth direction introduces the HV and VH components, and causes the shift of the 
location of the maximum response away from the orientation angle of 90°. The amount of shift correlates 
with the slope of the tilted surface. This phenomenon has been utilized to develop slope measurements of 
mountainous desert terrain and forested areas (Schuler, Lee and De Grandi, 1995). For topographic 
measurements, the inhomogeneities of forest and deserted areas often introduce various dominant 
scattering mechanisms. Algorithms for the slope measurement have to be carefully devised to account for 
this variation. For ocean surfaces, however, the dominant scattering mechanism for incident angles 
between 25° to 65° has been adequately represented by composite surface, or tilted-Bragg models for the 
backscatterer (Valenzuela, 1968). Consequently, the estimation of ocean surface slope is expected to 
produce more consistent results. 

This new slope measurement technique is applied to the P-, L-, C-band polarimetric SAR images 
of a current front with converging flows. The orientation angle displacement in the co-polarization 
signature is computed for each pixel. The results show a sudden 1.3° slope change at the front. All three 
bands show a similar trend. The effect of velocity bunching on this measurement technique will be 
investigated. The measurement error and the limitations of this technique will also be discussed. 

2. IN-SITU MEASUREMENTS 
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During the 1990 Gulf Stream Experiment, the AIRSAR imaged the north edge of the Gulf stream 
supersite near the coast of Virginia (36.5°N, 72.0 °W) ( Valenzuela , 1991). Simultaneous shipboard 
measurements of currents, temperature, salinity, etc. were also made. Two relatively bright lines in the 
SAR imagery delineated a filament of cool, fresh water that was entrained by the Gulf Stream as it passed 
near the continental shelf. In-situ measurements made from the R/V Cape Henlopen are summarized in 
Fig.l. A southerly bright line corresponds to a surface front having relatively weak current shear, but a 
convergence velocity (measured at 7-m depth) of about 20 cm/s. Nearly coincident interferometric SAR 
measurements confirm the weak shear and show a surface convergence of 20-25 cm/s. Surface 
convergence was also suggested by the narrow width of the frontal temperature change, surface 
accumulation of Sargassum, and an adjacent band of breaking waves. A more northerly bright line, which 
was not well-sampled, appeared to coincide with a zone of stronger current shear between the filament and 
the slope water to the north. Some convergence may have existed at this shear front as well. A second 
filament lying farther inside the Gulf Stream also showed a strong linear feature and also had strong 
surface convergence. Ocean waves were driven by southwesterly 8 m/s winds. Wave-current interactions 
at these fronts produced narrow bands of steep and breaking waves manifesting higher radar returns in 
polarimetric SAR images. Detailed discussion of the environmental setting and in-situ measurements is 
given by Marmorino et al (1994). 


(a) PLAN VIEW 



converging current front 
(b) VERTICAL SECTION 



Fig. 1. Schematic summary of the observations. Plan view (a) shows the east-west surface front 
between Gulf Stream and low-salinity filament. Slope- Water lies farther to the north. Vertical 
section (b) indicates relative cross-front current motion and surface convergence. Water-current 
interaction produces narrow band of steep and breaking waves appearing in the SAR imagery. 
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3. POLARIMETRIC SAR MEASUREMENTS 

The AIRSAR imaged the north edge of the Gulf Stream with quad-polarizations (HH, HV, 

VH.and W) of P-Band (68 cm), L-Band (24 cm) and C-Band (5.7 cm) simultaneously. For each band, 
the four polarization channels are processed to form SAR images, and then converted to the Mueller 
matrix. For speckle reduction and data compression, four Mueller matrices of pixels in the azimuthal 
direction are averaged to form 4-look processing. The pixel spacing is 12.0 meters in azimuth and 6.7 
meters in slant range. The typical image covers an area of 12.4 km in azimuth and 9 km in slant range. 
The incidence angle varies across an image from 10 to 65 degrees. We also obtained two scenes using a 
newer JPL 16-look processor. The 16-look data provides higher image quality with better focusing and 
less speckle noise. The resolution remains the same as in the 4-look cases. Our analysis has concentrated 
on the 16-look data. 

Several flight paths of the AIRSAR crossed the convergent front several times, providing valuable 
data for analysis. A north-south pass was investigated (AIRSAR data: CM3597), and it was found that the 
cross-pol (HV) signature of the front are much higher than eo-pol responses (W and HH), and that P- 
Band HV has the highest signal-to-background ratio. To quantify the radar cross sections, a cut wa made 
across the convergent front at an incidence angle of 47°. Fifty samples are averaged to reduce speckle. 
Fig. 2 shows the radar cross section, a, for all three frequencies and for all three polarizations. Cross 
sections at W are the highest, 3-8 dB lower at HH, and 9-18 dB lower at HV. For all three bands, the 
HV signal-to-background ratio, da, for the convergent front is much higher (—2 dB) than those of HH 
and W (=s 0.5 dB), and P-Band HV has the highest da, followed by L-Band and C-band. The spatial 
width of the da peak is wider for P-band than for L- and C-band, but the peak value for P-band is less 
than for L-Band. These effects are caused by the longer SAR integration time for P-band. When the 
integration time is longer than the ocean wave coherence time, it causes smear in the azimuthal direction. 

It is interesting to note that P-band W has a very weak signature. The presence of strong signatures in 
HV polarization indicates that strong wave-current interaction is producing large surface tilts in the along- 
track direction. 


4. OCEAN SURFACE SLOPE MEASUREMENT 


The polarization signature was introduced by van Zyl et al (1987) to graphically represent the 
variation of scattering cross-section as a function of the ellipticity and orientation angles of the transmitted 
waves. It is a useful technique for describing the scattering mechanism of point and distributed targets. 
For idealized Bragg scattering, the scattering matrix can be represented by 


'Bragg 


a o' 

0 £> 


( 1 ) 


where a and b are real, and b>a (i.e., ff vv >ff hh ). The HV and VH polarization responses are zero. Using 
the polarimetric SAR data (CM3597), for incidence angles between 20° and 65°, co-polarization 
signatures of all three bands display the same characteristics as Bragg scattering. For incidence angles 
below 20°, specular scattering becomes dominant. 


For Bragg scattering, maximum co-polarization signatures occur at an orientation angle \p=90°, 
and ellipticity angle x=0° (i.e., at W-polarization) . However, a tilted ocean surface in the azimuth 
direction introduces HV and VH components, and causes a shift of the location of the maximum response 
in the co-polarization signature away from \{/-90°. The amount of shift from \j/= 90° correlates with the 
slope of the tilted surface. Intuitively, the maximum co-pol response occurs when the vertical antenna is 
perpendicular to the ocean surface in the azimuth direction. The surface tilt causes the shift of the 
maximum response away from the vertical position. The amount of shift can be visualized as the rotation 
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of the vertical-horizontal basis vectors to an angle where the vertical component is again perpendicular to 
the tilted surface. This phenomenon has been utilized by Schuler et a/. [13] for slope measurements of 
mountainous desert terrain, where the Bragg scattering mechanism is also dominant. 

This new slope measurement technique is applied to the same data of Fig. 2. to measure the slope 
change across the converging currents. The ideal Polarimetric SAR imaging alignment for this slope 



Fig. 2 To quantify the radar cross sections of the convergent front, a cut was made 
across the convergent front at an incidence angle of 47°. Fifty samples are averaged to 
reduce speckle. The radar cross section, a, for all three frequencies and for all three 
polarizations are shown. The azimuth pixel size is 12m, and the distance covered is 2.4 
km (200 pixels). Cross sections at W are highest, 3-8 dB lower at HH, and 9-18 dB 
lower at HV. For all three bands, the HV signal-to-background ratio, da, for the 
convergent front is much higher (= 2 dB) than those of HH and W (=* 0.5 dB), and P- 
Band HV has the highest da, followed by L-Band and C-band. 
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measurement is to have the front aligned in the range direction. The AIRSAR CM3597 image is one of 
the data that possesses this alignment. The P-Batid polarimetric data is used to compute the orientation 
angle of the maximum co-pol response for each pixel. The result for the lower right quarter of data is 
shown as an image in Fig. 3. This figure displays a sudden slope change encoded as a gray-scale change 
at the convergent front. The R/V Henlopen is crossing the front with its wake in a contrasting dark slope. 
To quantify the slope change at the front, a 50-line average below the ship (at incidence angle =47°) is 
plotted in Fig. 4A. It clearly shows that the surface slope is positive north of the front, negative or near 
zero south of the front, and a maximum at the front. The estimated slope at the front has an abrupt 
change from an average of +0.8° to -0.5° as shown in Fig. 7. In our experiments, we have found that 
the shift in orientation angle induced by the surface tilt depends on incidence angle and the cr^ to cr w ratio. 

Similar results were also obtained for L-band (Fig. 4B), except that it is noisier. C-Band with its 
shorter wavelength is more subject to local wave variations, the noise level is much too high. For this 
application, P-band polarimetric SAR is preferred. 


AIRSAR flight direction 



P-BAND 



Fig. 3 The polarimetric SAR measured surface slopes around the 
convergent front is shown as a gray-level image. The slope is computed 
for each pixel by the shift from 90° of the polarization orientation angle 
associated with the maximum co-polarization signature. The sudden 
change in slope at the front demonstrated the effectiveness of this 
algorithm. 
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Fig. 4 A 50-line average of the orientation angles at the dotted line of 
Fig. 3 (incidence angle=47°) is plotted. The orientation angle of 90° 
corresponds to 0° surface tilt. It clearly shows that the surface slope is 
positive north of the front, negative or near zero south of the front, and a 
maximum at the front. 


6. REMARKS: 

The procedure of finding the maximum polarization signature follows a specially devised scheme 
that accerates the computation and improves the accuracy. This procedure ( Schuler , Lee and De Grandi, 
1995) is based on the method of steepest ascent, and uses iterative methods to locate the maximum. It was 
found that this method can achieve the orientation angle accuracy of 0.01°. 

The polarization orientation angle shift for a given surface tilt is incident angle and HH to W 
ratio dependent. For a given incident angle and HH to VV ratio, a theoretically derived function relates 
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the shift of orientation angle to the surface tilts. This function is then applied to map the polarimetric 
SAR orientation angle shifts into surface tilts. 

7„ Conclusions 

The polarization and wavelength diversity of radar signatures associated with narrow Gulf Stream 
current fronts have been investigated using the P-, L-, C-Band airborne polarimetric SAR data. Using the 
polarization signatures, a new procedure has been formulated to estimate the change in mean slope across 
the front. NASA/JPL AIRSAR data is used for demonstration. 
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INTRODUCTION 

The NASA/JPL airborne SAR (AIRSAR) system operates in the fully polarimetric mode 
at P-, L- and C-band simultaneously or in the interferometric mode in both L- and C-band 
simultaneously. The system became operational in late 1987 and flew its first mission 
aboard a DC-8 aircraft operated by NASA's Ames Research Center in Mountain View, 
California. Since then, the AIRSAR has flown missions every year and acquired images 
in North, Central and South America, Europe and Australia. 

In this paper, we will briefly describe the instrument characteristics, the evolution of the 
various radar modes, the instrument performance and improvement in the knowledge of 
the positioning and attitude information of the radar. In addition, we will summarize the 
progress of the data processing effort especially in the interferometry processing. Finally, 
we will address the issue of processing and calibrating the cross-track interferometry 
(XTI) data. 


INSTRUMENT CHARACTERISTICS 

In AIRSAR, transmit polarization diversity is achieved by alternately transmitting the 
signals using horizontal or vertical polarizations. Receive polarization diversity is 
accomplished by measuring six channels of raw data simultaneously, both H and V 
polarizations at all three frequencies. The video data are digitized using 8-bit ADCs, 
providing a dynamic range in excess of 40 dB and, together with navigation data, stored 
on tape using high density digital recorders . The AIRSAR system also includes a real- 
time processor capable of processing any one of the 12 radar channels into a scrolling 
image. In addition to checking the health of the radar, the scrolling display is also used to 
ensure that the correct area has been imaged. Table 1 provides a summary of the 
AIRSAR system characteristics. AIRSAR can be operated in many different modes due 
to the complexity and flexibility of the instrument. The evolution of these radar modes is 
summarized in the following section. 
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Table 1. Summary of AIRSAR system characteristics. The parameters in ( ) apply to 40 
MHz chirp bandwidth configuration. 



P-band 

L-band 

C-band 

Chirp Bandwidth (MHz) 

20 (or 40) 1 

Chirp Center Freq. (MHz) 

438.75 

(427.5) 

1248.75 

(1237.5) 

5298.75 

(5287.5) 

Peak Transmit Power (dBm) 

62 

67 

60 

Antenna Polarization 

H/V dual microstrip 

Antenna Gain (dBi) 

14 

18 

24 

Azimuth Beamwidth (deg) 

19.0 

8.0 

2.5 

Elevation Beamwidth (deg) 

38.0 

44.0 

50.0 

Antenna Size (m) 

0.9 x 1.8 

0.5 x 1.6 

0.2 x 1.4 

ADC Sampling Rate (MHz) 

45 (90) 

Data Rate (MB/s) 

10 

NE SigmaO (dB) 

-45 

-45 

-35 1 

Nominal Altitude (m) 

8000 

Nominal Velocity (Knots) 

450 | 

PRF/Polarization Channel 

1 (programmable) x ground speed in Knots 

Slant Range Resolution (m) 

10(5) 

Azimuth Resolution (m) 

1 

Ground Range Swath (km) 

10-15 


RADAR MODES 

When AIRSAR flew its first scientific mission in 19S8, it was capable of imaging sites 
in P-, L-, and C-band simultaneously in polarimetric mode or L- and C-band along-track 
interferometric (ATI) mode. ATI mode was successfully used to image ocean currents 
and waves moving in the radar line-of-sight direction. Since then, more antennas and 
antenna switching networks have been installed to accommodate cross-track 
interferometric (XTI) and bistatic modes. XTI mode was successfully used to generate 
topographic maps of areas of interest whereas the bistatic mode was successfully used to 
collect data in conjunction with ERS-1 (CVV) and SIR-C (LVV and CVV). Figure 1 
shows the relative location of all the antennas currently available on the DC-8 and their 
polarization. Table 2 summarizes the evolution of these radar modes. 

As shown in Table 2 and Figure 1, the baseline of the C-band ATI mode was shortened 
significantly in 1991 by pairing up the newly added C-bt antenna with the C-SAR 
antenna in an effort to increase the sensitivity to shorter decorrelation time of ocean 
currents. Prior to 1995, the single frequency XTI mode (XTI1) was operated with one 
transmit antenna that provides the best possible SNR. The reason is that the transmit path 
via the C-tp antenna is 3 m shorter than that of the C-bt antenna, hence giving us better 
SNR. Since 1995, we have been experimenting with alternating the transmit antenna 
between the top and the bottom antennas. This effectively doubled the baseline and 
initial data analysis showed that the longer baseline produced DEMs (digital elevation 
models) with reduced RMS height error as expected. In addition, the newly added L- 
band XTI mode produced DEMs of slightly higher RMS height error due to shorter 
baseline length (scaled by wavelength) compared to those of C-band XTI mode, 
although much work remains to be done to calibrate the L-band XTI mode. 
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To produce accurate DEMs, we need to know the baseline precisely. To do this, we have 
also upgraded the Inertial Navigation System (INS) and the Global Positioning System 
(GPS) receiver in order to have more accurate knowledge of the location and attitude of 
the antennas. The upgrades are described in the next section. 
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C-SAR (H/V) 



Figure 1. View of relative location of the antennas (not to scale). 


Table 2. Summary of the available radar modes for AIRSAR. * Since 1994, P-band is 
allowed to transmit in the U.S only when the radar is in 20 MHz chirp bandwidth mode. 
If special clearance is obtained prior to the flight, P-band is then allowed to transmit in 40 
MHz chirp bandwidth mode also. 


Mode 

Date 

P-band* 

L-band 

C-band 



TX RX 

TX 

RX 

TX 

RX 

POLSAR fauad-Dol) 

1988 - present 

mm 


IiiAVliiMi 

ATI 

1988 - 1990 

P-SAR (H/V) 

L-fwd and 
L-SAR 

C-fwd and 
C-SAR 

ATI 

1991 - present 

P-SAR (H/V) 

L-fwd and 
L-SAR 

C-bt (V) and 
C-SAR 

XTI1 

1991 

P-SAR (H/V) 

L-SAR (H/V) 

C-bt 

C-bt/ 

C-tp 

XTI1 

1992 - 1995 

P-SAR (H/V) 

L-SAR (H/V) 

C-tp 

C-bt/ 

C-tp 

XTI1 - ping pong 

1995 

P-SAR (H/V) 

L-SAR (H/V) 

C-bt and C-tp 

XTI2 

1995 

P-SAR (H/V) 

L-tp 

L-tp / 
L-bt 

C-tp 

C-bt/ 

C-tp 


1995 

P-SAR (H/V) 

L-tp and L-bt 

C-bt and C-tp 
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NAVIGATION SYSTEM 

The original navigation system of AIRSAR consisted of a Honeywell INS with a ring 
laser gyro that determined the attitude of the aircraft and a Motorola Eagle 4-channel 
GPS receiver that provided the positioning information (latitude and longitude) of the 
aircraft. As technology advanced and our need for more accurate positioning and attitude 
information became more stringent, we purchased a new Motorola Six-Gun GPS receiver 
and a new Honeywell Integrated GPS and INS (IGI) in 1994. The Six-Gun GPS receiver 
has six channels and a much more stable clock compared to the old unit and provides 
positioning accuracy of 100 m using CA code. This receiver was integrated in the radar 
in 1994. The Honeywell IGI has a smaller and more sensitive ring laser gyro integrated 
with a GPS receiver capable of receiving the more accurate but restricted Precise 
Positioning Service (PPS) data. The specifications on this unit are: 0.02° heading 
accuracy, 0.01° roll and pitch accuracy, 0.03 m/s velocity accuracy per axis, and 16 m 
positioning accuracy with PPS. The IGI was installed on the DC-8 in 1994 but the data 
were recorded off-line and were not available in the radar header until the 1995 flight 
season. 

In addition, we have also experimented with differential GPS by using a Turbo Rogue 
GPS receiver on the aircraft in conjunction with another Turbo Rogue receiver on the 
ground. This experiment is usually supported by the GPS experts from another section at 
JPL and requires special post-processing to obtain positioning accuracy of better than 1 
m. 


DATA PROCESSING 

A variety of processors and processing techniques are utilized to process AIRSAR data to 
imagery. A real-time correlator is part of the AIRSAR radar flight equipment (the 
Aircraft Flight Correlator) and is used to produce low resolution (approx. 25 meter) two 
look survey imagery. The same on-board equipment is used to generate a slightly higher 
resolution (15 meter), 16 look image of a smaller area (12 km x 7 km) within 10 minutes 
of acquisition using the quick-look processor. These on-board processors are useful for 
assessing the general health of the radar and the success of data taking in real time. 

Final processing of selected portions of the data to high quality, fully calibrated image 
products happens in the weeks and months following a flight campaign. Currently, users 
may request images from two different operational processors, the synoptic processor and 
the frame processor. In the synoptic processor, the user specifies three data channels to 
be processed. About five minutes of raw data from each of the three selected channels 
are processed to 16 looks and amplitude-only image strips, covering about 40 km along 
track. In 40 MHz mode, the image strips would be 8 looks and 20 km long. These image 
strips cover about 9 km in the slant range direction for the 20 MHz mode and 4.5 km for 
the 40 MHz mode. 

In terms of frame processing, we currently support two processor versions: the AIRSAR 
processor (version 3.5x) and the new integrated processor (version 5.x), which is still 
under development (especially in XTI calibration). 1995 has been a transition year for 
the AIRSAR ground processing facility. The new integrated processor was developed 
mainly to process XTI data routinely since XTI mode has become increasingly popular. 
In order to do so, we needed a new processor that tracks and compensates for the motion 
of the aircraft since uncorrected motion translates into baseline error between the two 
antennas, which results in height error in the DEM. In addition to motion compensation, 
the new integrated processor is also capable of generating images with full range swath as 
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opposed to half range swath with the version 3.5 processor. The XTI processor still 
needs an accurate algorithm to determine absolute phase. In addition, better calibration 
is required to remove systematic height errors in DEM. 

As with the previous version, the integrated processor processes one minute of raw data 
of all available data channels into absolutely calibrated images in compressed Stokes 
matrix format that contains all the polarization information. If C-band cross-track 
interferometer data are available for the data take, the integrated processor will generate a 
digital elevation model and a local incidence angle map. By using the local incidence 
angle map, all output images will be geometrically and radiometrically corrected taking 
the topography into account and resampled to ground range with a 10 m by 10 m pixel 
spacing. The output images cover about 10-12 km in the range direction by about 10 km 
in the along-track direction for the 40 MHz mode, and about 20 km in the range direction 
by about 10 km in the along-track direction for the 20 MHz mode. Although the radar 
data rate allows us to image about 20 km in range swath for the 20 MHz mode, the 
increasing phase noise due to decreasing SNR as a function of incidence angle reduces 
the correlation between the two antenna channels. As a result, the RMS height error can 
be quite large in far swath due to poor SNR. 


DATA CALIBRATION 

The calibration of polarimetric data is well understood. Briefly, with the calibration tone 
in the receive chain and comer reflector verification, we are able to consistently produce 
polarimetric images with better than 3 dB absolute accuracy, better than 1.5 dB relative 
accuracy amongst the 3 radar frequencies, and better than 0.5 dB between the 
polarization channels. The relative phase calibration between the HH and W channels is 
better than 10°. 

The calibration of XTI data is much more challenging because various parameters, such 
as baseline vector, are involved in the XTI data processing. The absolute phase must be 

known in order to derive height information from the interferometric data without 2k 
ambiguity. The differential phase (between two channels) of the radar can be a function 
of system temperature. Therefore, we need to determine both absolute and differential 
phase for each data take. In addition, accurate knowledge of the baseline between the 
two antennas to a few millimeters is necessary to generate accurate DEMs. We have 
successfully used the comer reflector array at Rosamond Dry Lake to determine the 
baseline for C-band antennas and are currently working on the L-band antennas that we 
started operating in 1995. 


SUMMARY 

In this paper, we described the AIRSAR instrument characteristics, the evolution of the 
various radar modes, and improvement in the navigation system. In addition, we 
summarized the progress of the data processing effort and briefly addressed some of the 
challenges in calibrating the XTI data. We hope to resolve the phase calibration issues 
with the 1995 dual frequency XTI data in the near future so that we can provide users 
with DEMs at L- and C-band routinely. 
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1. INTRODUCTION 

The capabilities of using remote sensing data, and in particular multifrequency/multipolarization SAR data, like 
AIRSAR, for the retrieval of surface parameters, depend considerably on the specificity of each application. The 
potentials, and limitations, of SAR data in ecological investigations are well known. Because the chemistry is a 
major component in such studies and because of the almost lacking chemical information at the wavelengths of SAR 
data, the capabilities of using SAR-derived information in such studies are considerably limited. However, in the 
case of surface energy/water balance studies, the determination of the amount of water content, both in the soil and 
in the plants, is a major component in all modeling approaches. As the information about water content is present in 
the SAR signal, then the role of SAR data in studies where water content is to be determined becomes clearly 
predominant. 

Another situation where the role of SAR data becomes dominant over other remote sensing systems is the case 
of dense canopies. Because of the penetration capabilities of microwave data, which is especially superior as 
compared to optical data, information about the canopy as a whole and even the underlying soil is contained in the 
SAR data, while only the top canopy provides the information content in the case of optical data. In the case of 
relatively dense canopies, as has been demonstrated in this study, such different penetration capabilities provide very 
different results in terms of the derived total canopy water content, for instance. 

However, although all such capabilities are well known, unfortunately there are also well known limitations. 
Apart from calibration-related aspects (that we will not consider in this study), and apart from other intrinsic 
problems (like image noise, topographic corrections, etc.) which also significantly affect the derived results, we will 
concentrate on the problem of extracting information from the data. Even at this level, methods are still not fully 
well established, especially over vegetation-covered areas. 

In this paper, an algorithm is described which allows derivation of three fundamental parameters from SAR 
data: soil moisture, soil roughness and canopy water content, accounting for the effects of vegetation cover by using 
optical (Landsat) data as auxiliary. Capabilities and limitations of the data and algorithms are discussed, as well as 
possibilities to use these data in energy/water balance modeling studies. 

All the data used in this study were acquired as part of the Intensive Observation Period in June- July 1991 
(European Multisensor Aircraft Campaign-91), as part of the European Field Experiment in a Desertification- 
threatened Area (EFEDA), a European contribution to the global-change research sponsored by the IGBP program 
(Bolle et al., 1993). 

2. PARAMETERS WHICH ARE REQUIRED IN 
SURFACE ENERGY/WATER BALANCE STUDIES 

Although the actual parameters which are required in surface energy/water balance studies depend very much 
on the kind of modeling approach adopted in each case, it is a general agreement that accounting for the amount of 
water available, and changes in water content, both in the soil and in the plants, is always a major component, not 
only for water balance but also for the partitioning of available energy into latent and sensible heat flux components. 
However, the way in which each parameter enters into the model and the assumptions made by each model are 
always conditioning the so called 'sensitivity' to such model parameters. For this reason, accuracy requirements on 
the retrievals of each parameter can not be easily stated. 

The model used in previous studies (Moreno et al., 1994), which was actually a derivation from the Biosphere- 
Atmosphere Transfer Scheme (BATS) (Dickinson et al., 1993) with significant modifications and additions, used a 
total of about 70 parameters (many of them fixed to default values), from which about 15 are potentially derivable 
from remote sensing data, and about half of them directly from SAR data or by combination of SAR data with 
optical data. Such parameters are: (top-)soil moisture, soil roughness, canopy water content, Leaf Area Index, 
vegetation height (displacement height), Stem Area Index-Canopy 'roughness' (or canopy geometry parameters), soil 
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albedo as a function of soil moisture (see Fig. 4) r as well as other parameters indirectly derived from image- 
classification results. Other parameters (like surface temperature, cloudiness, etc.), are also used, but will not be 
discussed here. 

Because of the complexity of intervening effects, no definite limits can be put a priori for the accuracy 
requirements over each parameter, and then no fixed limits are put over the capability to provide soil moisture or 
other parameters from remote sensing data, partly because the use of this kind of data (spatial data) would also 
require re-parameterization in the models (the problem of handling spatially distributed data is another reason for the 
difficulties in deriving clear conclusions from ID-model sensitivity studies). As the final desired goals (errors of 
about 10 Wm -2 for the derived fluxes) are still far from the actual capabilities (including ground-based 
meteorological networks), all we can do is to try to achieve the maximum accuracy possible. Also, the present 
situation is that the models (including 3D models) use very poor surface parameterization because of the lack of any 
additional data, so that any information which can be provided from remote sensing systems (with all the involved 
limitations) would still be in any case very welcome. 

3. CAPABILITIES OF SAR DATA TO PROVIDE THE REQUIRED PARAMETERS 

The capabilities of SAR data to provide at least some of the parameters required by surface energy/water 
balance models, especially those related to water content, are well known, and actual use of these data has been 
made in previous field experiments (FIFE, EFEDA, HAPEX). The way in which this has been done is mainly 
through model-inversion techniques. However, the kind of model used (and we actually do not have an appropriate 
model for the behaviour of natural surfaces at the frequencies at which SAR data are acquired), and the kind of 
inversion technique used, become critical when retrieved numerical values are to be compared to ground 
measurements. As the SAR signal sensitivity to water content in the canopy, and also in the soil, has been in any 
case demonstrated (Engman, 1991), it is expected that this information will come out in the retrievals derived from 
SAR data, at least in terms of relative values. 

The use of full-polarimetry information for the retrieval of soil/vegetation parameters has been demonstrated to 
be an essential aspect as compared to the capabilities of single channel/single polarization systems (like ERS-1/2). 
Because of the ability of polarimetric information to separate different contributions, roughness/geometry effects 
and water content information can be decoupled, and then the corresponding values retrieved from the measured 
data. 

Apart from the classical parameters derived from polarimetric information, information about canopy roughness 
(mainly related to canopy height) can be also derived from interferometric information. However, because the 
physical meaning of canopy roughness (which is also a recurrent function of wind speed) cannot be easily related to 
such interferometric information, work in this field is still in its beginning. Phase information in standard 
polarimetric data can also be related to canopy geometry, but the link between such estimates and the required 
parameters (canopy height and displacement height, canopy 'roughness') is still at the level of empirical 
relationships. The penetration capabilities of SAR data become here a difficulty, because it is the canopy height 
which is expected to be given as input to the models. However, the actual canopy roughness is not simply related to 
vegetation height, and most probably the roughness information derived from the scattering mechanism in SAR data 
is more relevant for modeling purposes than the use of a fraction of the total canopy height as a roughness estimator, 
as is currently being done in the models. 

It is true that AIRSAR data alone still have limitations in this type of study because of the ambiguity in 
scattering mechanisms in (dry) bare soil surfaces and in surfaces with small vegetation amount. The combination of 
information derived from optical data (where the separability between soil and vegetation behaviour is more strict, 
see Fig. 3), with AIRSAR data, significantly improves the capabilities of the retrievals of the required information. 

Two main aspects are in any case to be taken carefully into account: data preprocessing aspects (becoming 
critical for any posterior analysis of the data), and the kind of information-extraction technique used to derive the 
required parameters from the data. Both will be discussed in the following sections. 

4. DATA PRE-PROCESSING 

In any study in which theoretical physically-based scattering models are to be inverted against measured data, 
the calibration of such measured data becomes critical. The reason is that the models are based on physical 
coefficients which are supposed to be of universal applicability. The way in which such parameters come into the 
model (usually through highly non-linear relationship) makes impossible the introduction of some kind of (linear) 
compensation for calibration deviations, so that irremediably any calibration error (including deviations in antenna 
gain pattern correction for varying altitude/topography and any other related radiometric correction) produces 
absolute errors in the retrievals by using the theoretical model. For the data used in this study, calibration was done 
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according to standard techniques at JPL. A total of four comer reflectors were deployed in the study area during the 
AIRSAR overflight, and their response was used to check calibration in the data and to perform the proper 
corrections when needed. 

The synthesized images were geometrically rectified to ground-range projection by using available ephemeris 
data for the DC-8 aircraft navigation, while at the same time azimuth/range pixel sizes were compensated for to 
make square the resulting pixels through cubic convolution. Finally, the images were resampled by using about 50 
ground control points in the common 10 m UTM grid, covering the 10x10 km 2 area, to which all remote sensing 
data were co-registered into the image database to make possible the use of multisensor/multitemporal studies. As 
in the resulting geometrically corrected image, pixel numbering loses the information of the incidence angle, a new 
image with the corresponding incidence angle for each new UTM pixel was also produced to facilitate additional 
processing. Although accurate Digital Elevation Models are available for the area, specially developed as part of the 
EFEDA experiment, because the study area considered here is completely flat (maximum height differences of less 
than 20 m and almost constant slope over the full area), no topographic corrections have been applied. 


5. MODEL INVERSION TECHNIQUES APPLIED TO MULTI-PARAMETER 
RETRIEVAL FROM SAR DATA 

The first approach that we considered in the derivation of surface parameters from SAR data was merely 
empirical, based on correlation approaches between some "indices" (band and/or polarization combinations) and 
measured surface parameters (LAI, biomass, canopy water content, etc). Different band ratios (mainly L/C and P/L) 
and polarization ratios (HH/W, HH/HV, W/HV) have been considered. Although not intended to be used for the 
actual retrieval of surface parameters from the data, the establishment of empirical relationships allows a preliminary 
estimation of the capabilities of the data to account for the observed variability (and the determination of error 
bounds which can be expected and used in the fitting of the 'merit function’ in the numerical inversion procedure) as 
well as the derivation of simple relationships to be used in the initialization of the model parameters in the iterative 
inversion procedure used later for parameters’ retrieval. Although the advantages of empirical relationships are well 
known in terms of speed in calculations and avoidance of convergence and other numerical problems, the use of 
empirical relationships is absolutely limited by how well one can extrapolate from the results and the generality of 
the algorithms. Then, even if the use of empirical relationships could provide an optimum fit in our case where we 
have ground measurements of all surface parameters, the use of such empirical methods has been avoided and they 
have been restricted to be used as auxiliary elements in the numerical model-inversion procedure. 

The central part of this study is the development of a model-inversion technique to extract the required 
information from SAR data. The details of the method are given elsewhere (Saatchi et al., 1993; Moreno et al., 
1994; Moreno, 1995; Moreno and Saatchi, 1996), and only the main aspects will be described here. 

For the contribution of bare soil, several models have been considered in this study. A model developed by 
Dubois et al. (1994) was initially used (Saatchi et al., 1993). Another semi-empirical model developed by Oh et al. 
(1992) has been also used. The present implementation is an updated version of the Oh et al. algorithm after more 
recent improvements introduced by the same authors. The results obtained for bare soil moisture are indeed in 
agreement with the results obtained by other authors (Oevelen et al., 1995) by using the Integral Equation Model 
over the same dataset. As a major effort is put on the derivation of canopy parameters, because the model is 
intended to be used over agricultural areas with significant vegetation cover for some fields, the parameterization of 
the soil has to be kept to a minimum in order to make the model actually invertible. 

For the derivation of canopy parameters, the model has been implemented with a layer of randomly distributed 
scattering elements over the underlying bare soil. Such a layer represents the vegetation contribution. Three types 
of scattering mechanism are then considered: volume scattering, surface-volume scattering and soil surface 
scattering. The direct scattering from the soil is also attenuated (twice) because of the presence of the canopy, and 
such attenuation must be also included in the model. The surface-volume scattering terms must be introduced 
because of the type of vegetation to be considered in this case, giving significant contributions only for the co- 
polarized terms. The backscattering coefficients for the canopy are obtained by using the distorted Born 
approximation (Lang and Sidhu, 1983). Finally, the model needs to account for the amount of effects which are due 
to the soil and those which are due to vegetation in the case where no dense vegetation is considered but sparse 
vegetation is (and is still assumed randomly distributed at the scale of a pixel in order to avoid problems in the 
modeling of very clumpy structures). 

An important aspect to be pointed out is the necessity of working with the original channels separately in order 
to minimize the problems originated by the presence of noise and in order to take full information from the data, 
without reductions in the dimensions of the original information. Most of the empirical approaches work with 
channel ratios. These channel ratios are supposed to somehow compensate for deviations in calibration and 
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secondary effects in the signal, but this is not absolutely true, and also these channel ratios actually enhance the 
noise as compared to the noise present in each channel separately. The potential advantages of using channel ratios 
do not compensate for the problems that the use of channel ratios introduces in the inversion procedure. 

The method used for numerical inversion of the scattering model is the downhill simplex method, with two 
limiting conditions: maximum error and maximum number of iterations allowed. Although less robust (and 
especially more time consuming) than other techniques, this method has proven to be more resistant to noises and 
inadequacies in the model to fit the data, providing always a set of solution parameters for each pixel after avoidance 
of divergences. 

6. RESULTS 

Although the data used in this study correspond to only one specific situation (which prevents us from deriving 
general conclusions), and also the data is old as compared to the new capabilities added to the AIRSAR system in 
the last few years (especially in terms of data calibration), the conclusions derived from this study are similar to 
those previously derived in other field experiments: techniques are promising but still not fully ready for operational 
application, at least for the purpose of energy/water balance monitoring. 

In the case of soil moisture, because of the very high sensitivity of energy/water balance models to the 
initialization of the soil moisture profile prior to temporal evolution calculations, the values retrieved from SAR data 
are not very useful for this purpose. However, those values are in good agreement with ground observations. This is 
especially significant over the study area used in this case, because it is a very dry area (especially in summer, when 
the experiment was carried out) and it still seems that the SAR signal is sensitive enough to the top-soil moisture 
content, provided that the effects of surface roughness variations are properly accounted for. 

In the retrieval methods based on single-channel multiple-polarization techniques we have observed that in 
many cases the retrieved values of soil moisture at C band are considerably higher that those at L band, and this is 
independent of the kind of model inversion used. Because this is in contradiction with the expected behaviour 
(taking also into account the measured profiles given in Table 1), we have developed a double-channel multiple- 
polarization inversion technique (which is also actually necessary in order to account for vegetation effects in the 
retrievals). In the multiple-channel approach, the results are more consistent. Also, the resulting surface roughness 
terms are more consistent when dual channel methods are applied, especially for the correlation length. However, 
the use of dual-channel approaches requires that the model be applicable to both frequencies, and this is not true for 
some roughness conditions. In any case, the uncertainty in the retrievals of soil moisture (at least over our study 
area) is about 20-50% for the top-soil moisture content. 

As part of the data collection for the field experiment, the soil group generated a detailed soil map of each pilot 
area, including soil type, soil texture, soil depth and other information. Soil density and hydrological properties 
were measured both in the field and in laboratory conditions, and special experiments were carried out to test the 
importance of spatial variability in such soil properties, even at different spatial scales, from the field level up to a 
network of 10x10 km 2 (Bolle et al., 1993; Droogers et al., 1993; Ogink-Hendriks et al., 1995). However, if all this 
information is used to calculate the dielectric constant of the soil for comparison with the retrievals derived by 
model inversion from SAR data, we are facing the problem of empirical relationships between soil moisture and 
dielectric constant (Hallikainen et al., 1985). Because of the large uncertainty in such relationships over varying 
natural conditions in the field (also coupled to high variability in surface roughness as modeled by statistical 
estimators based on rms heights and correlation lengths), absolute values of soil moisture must be regarded with 
some high probability of 'systematic' errors. However, relative values are in any case comparable, but for that 
purpose an independent account of the variability in intrinsic soil properties (texture, density, etc.) is needed, and 
this is something which is rarely available. Even when such information is available for pilot areas (as in the case of 
the EFEDA experiment), the way in which such information can be used in conjunction with SAR data remains 
unclear because of the problems of spatial scaling and scale compatibility. 

In the case of canopy water content (see Fig. 2), the results are very sensitive to the geometric characterization 
of the canopy and the soil. The retrieved values (see for instance Fig. 2) are quite reasonable, especially taking into 
account that no ground data are used for training the model but only theoretical scattering considerations. However, 
the model was developed for the case of corn canopies in particular, and, in the present version of the algorithms, 
extrapolation to a full image requires some pre-classification of the scene in order to account for scattering 
mechanisms in a different way over different vegetation types. Work in progress is trying to eliminate this 
dependence by introducing additional parameters accounting for canopy geometry effects, but it seems that such 
geometrical effects can only be accounted for by means of multiple data taken with different incidence angles. 
Otherwise, separation of geometric effects from actual canopy water content will never be possible for absolute- 
value retrievals. 
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The modeling of soil roughness underlying the canopy has been demonstrated in this study to have also a major 
importance. An underestimation of soil roughness results in an overestimation of canopy water content. 

The first trials to get information from SAR data alone resulted in difficulties in the interpretation of the 
retrieved values in the case of partially covered pixels. Because the soil/vegetation algorithms used in combination 
worked only over bare soil areas or over dense (homogeneous) vegetated areas, the cases where separation between 
both extreme cases was not quite obvious give wrong contributions for soil moisture and/or canopy water content. 
The introduction of Landsat TM data (see Fig. 1) as auxiliary information, and the use of just the fractional 
vegetation cover from optical data (instead of other parameters which could be also derived from Landsat data) give 
as a result a significant increase in the capabilities of the application of the algorithm (in terms of reduction of the 
number of iterations needed and avoidance of cases with no convergence as always forcing a linear solution). 
However, because of the difficulties in modelling partially vegetation-covered soils, the results are still questionable 
in such cases. 


7. CONCLUSIONS 

A key point in the results is the necessity of some spatial homogenization of the original data prior to any 
information-extraction technique being applied, especially in the cases where model-inversion techniques are used. 
Multilooking techniques (spatial average) have been demonstrated to be not enough for the purpose of deriving 
consistent moisture fields from the data. It is necessary to reduce the spatial resolution to about 100 m to get 
consistent fields over homogeneous areas, but, as the spatial resolution decreases, more difficulties are added in the 
inversion technique to get convergence over heterogeneous areas. Actually, as spatial resolution is decreased, the 
noise level is reduced, but some information is lost. However, the main impact of reducing spatial resolution is 
increasing the within-pixel heterogeneity and decreasing between-pixel variance. The result is that over 
heterogeneous pixels the model cannot be inverted because no convergence is possible or because the retrieved 
values are out of range for many of the resulting heterogeneous mixtures (lack of physical meaning for the model). 

In order to get reliable results it is necessary to work with the highest possible spatial resolution but allow a two- 
way accounting for inter-pixel variability which is just due to noise: spatial filtering and multiresolution inversion 
techniques. The combination of both techniques is the only realistic way to handle the problem of spatial variability, 
especially when one of the objectives is just to analyze the problem of spatial variability in the derived surface 
values. 

The capabilities to derive soil moisture values which can be realistically used in surface energy/water balance 
studies seem to be very limited, especially over dry areas, where no sensitivity to deep moisture content is present in 
the SAR signal but where vegetation can take water from very deep levels and still produce a considerable 
contribution of latent heat flux in surface energy partitioning. According to our results (see Tables 2 and 3), an 
uncertainty in soil moisture values between 20-50% is all we can get in the case of relatively dry areas. Note 
however that in situ measurements in controlled conditions also give uncertainties between 16% and 40% (Table 1), 
in the same range as the variability in the different values derived from SAR data. Even when top-soil moisture can 
be detected by SAR data, the use of these data in energy/water balance monitoring is still very limited in vegetated 
areas, because of the predominant role of root-zone moisture, as well as in the case of bare soil, where other 
mechanisms are determining more strongly the dynamics of water in the soil. Whatever these limitations, the SAR- 
derived top-soil moisture field can still be very useful to properly account for the variability of soil albedo as a 
function of soil moisture. In the case of soils with low albedo, variations in soil albedo due to changes in soil 
moisture can be up to 50% (Fig. 4). As the soil albedo has a major effect on energy balance (actually controlling the 
amount of energy which is available at the surface), any improvement in surface albedo retrievals (including 
temporal variability due to soil moisture changes) would have a dramatic impact in energy balance monitoring 
through modeling techniques, provided that an observation system, stable enough, could provide routine updates of 
systematic top-soil moisture changes. 

Concerning canopy water content, the results obtained from AIRSAR data seem to overestimate the measured 
values of the ground according to the model used in this case. Although there are several uncertainties in the model 
that could explain such differences (especially those related to canopy geometry factors), it seems that the critical 
aspect in the model is the way in which soil scattering is treated. The overestimation in canopy water content is then 
mainly due to an underestimation of the roughness of underlying soil. It is interesting, on the other hand, to compare 
the results obtained from AIRSAR data to those obtained from optical data (Landsat TM) (see Fig. 2). In the case of 
Landsat TM, a full model inversion technique is also used. The model takes the reflectance values measured in the 
six (thermal channel 6 is excluded) Landsat channels (after radiometric calibration and atmospheric correction of the 
data). Then, an inversion technique is applied to fit the 6 measurements to give 6 surface parameters. One of the 
surface parameters used, and then retrieved by the algorithm, is the canopy water content. The values of canopy 
water content retrieved from Landsat data represent a strong underestimation of the values measured on the ground. 
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The reason for that is the sensitivity of optical data only to the top of the canopy (top leaves), so that a high 
underestimation is always expected from such optical data. However, fractional vegetation cover values derived 
from optical data are essential to guarantee proper interpretation of SAR data over partially covered areas or where 
confusion between soil and vegetated areas can be present in SAR data, as has been demonstrated in this study. 
Then, optical-microwave synergy seems to be the only way to overcome the limitations of both optical and 
microwave data in energy/water balance studies. 
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Table 1. Measurements of soil moisture over a reference bare soil field used for multisensor calibrations/intercomparisons 
(located just at the center in Fig. lc), during the two AIRSAR overflights over the study area. Numbers correspond, respectively, 
to the mean (of the N values available), standard deviation and relative error (%). The N measurements correspond to a spatial 
grid over the field. 


Date 

Soil moisture 

Soil 

roughness 



TDR (N=29) 

Volumetric: 

a (cm) 

X (cm) 

19 June 91 (a) 



crust : 

2.911.2 (40%) 



M 

5 cm : 

3.5 ± 0.8 (24%) 

0-5 cm : 

4.910.9(19%) 

1.41 

4.00 

ii 

10 cm : 

7.2 ± 1.2(17%) 

5-10 cm : 

10.911.8(16%) 



n 

30 cm : 

16.013.5 (22%) 







Gravimetric: 

(N=24) 




14 July 91 (b) 

0-5 cm : 

4.0 + 0.7(17%) 





(a) Droogers et al., 1993; (b) Saatchi et al., 1993. 


Table 2. Derived soil moisture values and roughness parameters from AIRSAR data for the same reference bare soil where 
simultaneous ground measurements are shown in Table 1. 


Date 

Soil moisture 

Soil roughness 


L band 

a (cm) 

X (cm) 

19 June 91 

3.9 

0.19 

3.6 

14 July 91 

7.1 

0.27 

3.1 


Table 3. Comparisons of retrieval of soil moisture by using the same AIRSAR data but applying three different models for the 
same reference bare soil area where simultaneous ground measurements are shown in Table 1. Data shown correspond to model 
inversion for L band (HH-VV and HH-VV-HV, depending on the model used in each case). 


Date 

Soil moisture 


Ohetal., 1992 

Fung et al., 1992 

Dubois et al., 1994 

19 June 91 

11.2 

5.6 

6.8 


63 



















(c) Landsat T 
ill dataset. T1 




BARRAX-BM1 field / 14 July 1991 



Fig. 2 Comparison of the retrievals of canopy water content from almost simultaneous microwave (AIRSAR) and 
optical (Landsat TM) data. SAR data slightly overestimate the total canopy water content in this case, but Landsat 
data give a very low value as compared to ground measurements. 



Fig. 3 Backscattering coefficients (C-HV) as measured by AIRSAR plotted against NDVI values derived from 
simultaneous Landsat TM data, for all the pilot fields used in the EFEDA’91 experiment (MAC-Europe campaign). 
Optical data allows a better separation of bare soil and vegetated surfaces for interpretation of SAR data. 



Surface moisture 

Fig. 4 Bare soil albedo (spectral and angular integration) as a function of top-soil moisture, according to the 
parameterization used in the Biosphere-Atmosphere Transfer Scheme (BATS) surface energy/water balance model 
(numbers on labels indicate the soil albedo corresponding to saturated conditions for each curve). 
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AN AUTOMATED MAPPING PROCESSOR USING 
C-BAND INTERFEROMETRIC SAR DATA 


E. Rodriguez, T.R. Michel, J.M. Martin, B. Houshmand 
Jet Propulsion Laboratory 
California Institute of Technology 

Abstract: We present the description of a processor which has been implemented to 
generate map products starting from C-band interferometric data. The first stage of the 
processor consists of the conventional interferometric SAR processing producing a Digital 
Elevation Model (DEMs) and a SAR brightness image in sensor coordinates. In the second 
stage of processing, a land use classification map is obtained by using the DEM, 
brightness, and interferometric correlation layers. Auxiliary layers which include a drainage 
layer, a height gradient layer, a height error layer, an estimated penetration layer, and a 
shaded relief layer are also computed. In the final step, all UTM collocated layers are 
combined in a GIS system which allows for both hard copy map products and for digital 
applications. 


1. INTRODUCTION 

Conventional generation of map products using stereo imaging techniques requires human 
intervention, making the elaboration of wide maps covering very extensive areas a daunting 
proposition. Recently, Interferometric SAR (IFSAR) techniques have shown the capability 
of producing digital elevation models (DEM) in near real time, and without operator 
intervention. Given the promise shown by these techniques, it is a natural extension to try 
to implement a processing system capable of generating digital map products automatically 
in near real time. The purpose of this paper is to report on such a prototype system which 
has been implemented at JPL for the processing of AIRSAR TOPSAR C-band 
interferometric data. 

The goals in developing the automated mapping station described here were the following: 

1) Automated product generation without operator intervention or editing. 

2) Near real time processing capabilities. 

3) Land classification from interferometric data alone. 

4) Digital products in a GIS usable format. 

5) Ability to produce hard copy output of standard map products. 

The following section describes the JPL AUTOMAP processor which has been 
implemented to achieve these goals. 

2. PROCESSOR IMPLEMENTATION 

The AUTOMAP processor takes raw interferometric data (e.g., the data generated by the 
JPL TOPSAR radar), and produces two types of products: 1) Geolocated topographic 
maps; and, 2) A multi-layer GIS database which can be examined on a work station. The 
AUTOMAP processor can be divided functionally into four (see Figure 1): 
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• The IFPROC Processor: Converts raw interferometric data into heights, radar brightness, 
and correlation files. 

® GIS Data Layer Generation: Takes the outputs from IFPROC and computes the layers 
which will be incorporated into the GIS database. 

• Geolocated GIS Quad: Given the desired map location and size, creates a geolocated 
database incorporating the desired data layers. 

• Map Product Generation: Makes map products suitable for printing. 

The IFPROC processor has been documented elsewhere, and this manual will restrict itself 
to a brief overview of how it is integrated to the rest of the AUTOMAP processor. 

The IFPROC module of the processor always outputs four data files which correspond to 
the height, brightness, correlation, and header data. Since user requirements for the output 
map products are variable, the rest of the AUTOMAP processor is designed for flexibility 
in both the output products and in the contents of the GIS database. The user can specify 
the desired products by using a menu, which can be filled interactively or preset prior to 
run time. 

The processor is currently implemented in a 16-processor Silicon Graphics Power 
Challenge computer and a DEC Alpha work station which is used as a front-end for 
displaying the GIS data and printing output products. The computer intensive code is run in 
the Silicon Graphics computer and is written in C++, C, and FORTRAN. The 
commercially available package PCI is used to geolocate the data layers, do some 
processing steps, and display the data. PCI processing steps are run automatically using 
the PCI "easi" script language. 

Terrain Classification 

The terrain classification module takes in the height, brightness, and correlation data 
produced by IFPROC, computes a series of data feature layers, and performs terrain 
classification using a mixture of Bayes classification and knowledge based algorithms. A 
flow diagram is presented in Figure 2. 

The following is a brief description of each of the programs in the terrain classification 
module: 

• featureVector: This program takes the three files generated by IFPROC and generates a 
set of user defined features to be used by the Bayes classifier. The features currently 
implemented are: 

• The radar brightness corrected for incidence angle, antenna pattern, and range 
effects. 

® The measured correlation between the two EFSAR channels. 

• The calculated "penetration" (correlation corrected for 

thermal noise and tilt effects). To obtain this feature, define the volumetric 
decorrelation, Yz> as 


* YgYn 

where Y is the measured correlation; Yn is the decorrelation due to noise, defined as 
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Yn = 


1 


1+ SNR 1 

where SNR is the signal-to-noise ratio; and, Yg is the decorrelation due to the size 
of the resolution cell and the local surface slope. The "penetration", Gz , which is a 
measure of the standard deviation of scatterers in the vertical direction, is defined as 

Oz = :fVl “ Yz 
Kz 

• The height rms inside a window whose size is determined by the user. 

• The slope rms inside a window whose size is determined by the user. 

• The brightness rms inside a window whose size is determined by the user. 

• The rms of the detrended slope. 

• The rms of the detrended height. 

In addition, featureVector calculates the surface slopes in the north and east directions. 

• bayesCIass: This module implements a standard Bayes classifier assuming Gaussian 
probability density functions and equal a priori probabilities. The inputs consist of the 
features calculated by featureVector and a file which contains the training information 
required by the Bayes Classifier. The philosophy behind this approach has been to 
decompose the true feature pdfs using Gaussian basis functions, and each basis function is 
treated as a separate class in the classification. For this version of the processor, a total of 
four land use classes have been chosen: water, trees, urban, non-treed/non-urban terrain. 

To go from the many (typically on the order of 20) classes used in the Bayes classification 
to these four classes, a projection method is used, as described below. Training was 
accomplished by using sites for which ground truth had been manually collected. The 
training set spanned a variety of sites including areas around the San Francisco Bay area, 
Crater Lake, Washington, and the Los Angeles County area. 

• collapseClass: As discussed above, the number of classes generated by bayesCIass 
is greater than the four classes ultimately used in the terrain classification layer. 
collapseClass uses a simple mapping to collapse the classification layer into four classes. 

• cleanClass: This program uses contextual information to clean up the water layer in the 
classification. The information used includes surface slope and area covered by connected 
components. 

• sieve: In general, the classification will have a greater resolution than is desired in the 
output map product. The sieve module decomposes the classification map into connected 
components and iteratively merges smaller ones into neighboring larger ones, until every 
component in the classification map has a minimum area. This procedure avoids the 
blockiness associated with simple averaging algorithms and preserves edges. 

Feature Calculation 

The AUTOMAP processor currently implements only two automatically derived feature 
layers: the mountain/peak and drainage layers. 

• Mountain/Peak Layer: This layer is generated by two programs, mountain and 
peak, mountain computes a mountain mask based on slope and covered area criteria. 
This mask, together with the elevation, is used by peak to find ail the major mountain peaks 
contained inside the mountain mask. 
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9 Drainage Layer: Due to its large memory requirements, the drainage layer is calculated 
at the map quad level, rather than for the strip map data. The drainage network is calculated 
by a module of the TOPAZ topography analysis package, version 1.10, produced by Martz 
and Garbrecht for the US Department of Agriculture. The complete TOPAZ package can 
measure topographic properties, define surface drainage, delimit watershed boundaries, 
quantify the drainage network, and parametrize subcatchments and overland flow paths, 
but we use only a small subset of its capabilities. 

To obtain a drainage network, the DEDNM module first preprocesses the DEM to remove 
localized depressions by filling them in up to the lowest outflow point, and lowering the 
DEM levels of localized impoundments (dams) across drainage paths to allow the flow to 
continue. Hat regions either in the DEM or produced by depression fill-in are then 
modified by adding small height variations consistent with the surrounding terrain and just 
large enough to define a flow path across the regions. The D8 (deterministic eight- 
neighbor) method is used to calculate the steepest descent direction for each pixel. Then, 
starting at each DEM cell in turn, the flow vectors are followed from cell to cell until the 
edge of the DEM is encountered, and the upstream area value for each cell along the path is 
incremented. After all DEM cells have been used as starting points, the upstream area value 
for each cell is the number of DEM cells that drain into that cell. Finally, to define the 
drainage network, only cells that have an upstream area greater than a threshold value (the 
Critical Source Area, CS A) are marked, and paths shorter than a min im um length 
(Minimum Source Channel Length, MSCL) are pruned. The CSA and MSCL can be 
spatially varying to account for varying hydrologic controls, such as soil type or 
vegetation, but we use only a single value for these quantities for the drainage network 
calculation. 


The DEDNM module has been slightly modified for the AUTOMAP purpose. The original 
program did not have an exit path after the drainage network was calculated, so one was 
added as an option. The DEM input routine was also modified to read a binary DEM file, 
and to accomodate data values less than lm. Finally, the output is reported as a vector 
layer, rather than a byte image. 

Auxiliary Layers 


In addition to the previous layers, AUTOMAP computes a series of useful auxiliary layers 
which may be incorporated in the GIS data bas and output map products: 


* Height Error Layer: A layer containing the random DTE height error, as derived from 
the measured correlation function, is also calculated. The height error is calculated by first 
calculating the estimated phase noise, 0$, using the Cramer-Rao bound 


Gg> “ 1 


vr 


f2N y 

where N is the number of looks. The height error is then estimated using the formula 

_ A r (sin 0 + cos 0tan x) 

2tzB 


where B is the component of the baseline perpendicular to the look direction, 0 is the local 
incidence angle, r is the range to the scatterer location, x is the surface slope in the azimuth 
direction, and A is the wavelength. 
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• Shaded Relief: A shaded relief image is calculated from the slope layers. The 

•* /S /\ 

brightness value is set by forming the dot product, n*r, where n is the local surface 
normal, and r is the unit vector in the look direction. 

• Brightness Modulated Classification RGB Layers: Red, green, and blue image 
planes are generated using the brightness to set the image value and die classification to set 
the image hue. 

• SigmaO Layer: The radiometrically corrected brightness feature is used as an auxiliary 
layer. 

• Penetration Layer: The penetration feature is used as an auxiliary layer. 

Geolocation and GIS Database Generation 

The geolocation is performed using nearest neighbor interpolation going from sensor 
coordinates to UTM. The GIS database generation is performed by using PCI, a 
commercially available GIS program. Given the user selected map size and location and a 
set of options, this module will resample the strip map layers generated in the previous 
module to a UTM map using the WGS 84 ellipsoid. It will, in addition, calculate a contour 
level vector layer. The output GIS data base will be stored in PCI format. 

Map Product Generation 

The map product generation module is performed by using PCI hardcopy map making 
utilities to generate the map data product selected by the user. The output is in a standard 
graphic format which can be output by a color printer or viewed electronically. 

3. SAMPLE RESULTS 

Figure 3 presents a sample output of the AUTOMAP processor for data over the San 
Francisco Area. The training set for this classification was obtained by selecting small 
representative samples for each of the four classes from the data itself. In the image, 

Golden Gate Park and the Presidio region, which are known to be treed, are clearly 
recognizable. The urban areas agree with USGS 7.5' quads. The Bay and Pacific Ocean 
are also clearly identified. A detailed comparison against ground truth confirms many of the 
detailed features observed in the classification map, although there are areas, such as in the 
hills above Golden Gate Park where the extent of tree coverage is over-estimated. The 
contour lines shown in the map are obtained from the TOPS AR IFS AR derived DEM. 

4. CONCLUSIONS 

We have presented a description of the JPL AUTOMAP processor for generation of map 
products and auxiliary layers using C-band IFS AR data done. Preliminary results indicate 
that the techniques presented here represent an automated alternative to traditional stereo 
photography for generating finished map products. Furthermore, the current processing 
speed of the AUTOMAP processor is 3 squared kilometers/per minute using a 16 processor 
Silicon Graphics Power Challenge. This efficiency is much greater than that which can be 
achieved using conventional optical techniques. 
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GIS LAYER GENERATION 



Figure 1: Flow diagram for the JPL AUTOMAP Processor 













4 ^ 



Figure 2: Flow diagram for the terrain classification module of the JPL 
AUTOMAP processor. 








Figure 3: Sample product of the JPL automap processor showing IFSAR generated 
height contour levels, terrain classification (image hue), and radar brightness (image value) 
simultaneously. 



San Francisco, CA — Automated Classification Map 
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First Results of the TOPSAR C-Band / L-Band Interferometer: 
Calibration and Differential Penetration 


Paul A. Rosen and Scott Hensley 


Jet Propulsion Laboratory 
California Institute of Technology 
4800 Oak Grove Drive 
Pasadena, CA 91109 


SUMMARY 

The NASA/JPL TOPSAR instrument recently was extended from a single- wavelength C-band 
(5.6 cm- A) dual aperture synthetic aperture radar interferometer to include a second wavelength 
at L-band (24 cm). Adding the second wavelength invites comparison of wavelength-diverse effects 
in topographic mapping of surfaces, with the principal goal of understanding the penetration of 
the radar signals in vegetation canopies, and determining the inferred topographic height. A first 
analysis of these data was conducted at two sites. Elkhorn Slough near Monterey, California 
presented flat, vegetation free terrain required for calibrating the radar interferometric parameters. 
A second site stretching from San Jose to Santa Cruz, CA, which is heavily vegetated, provided 
the first test case for wavelength diverse penetration studies. Preliminary results show that: a) the 
interferometer calibration determined at Elkhorn Slough is extendable to Laurel Quad and gives 
confidence in the C- and L-band height measurements; b) Clear differences are observed between 
the C- and L-band heights associated with vegetation, with C-band-derived topographic heights 
generally higher than those from L-band. The noise level in the L-band interferometer is presently 
the limiting factor in penetration studies. 

Introduction 

The NASA/JPL TOPSAR instrument, maintained and operated by the Jet Propulsion Laboratory 
AIRSAR group (Yunjin Kim, supervisor) was extended from a C-band-only topographic mapping 
SAR interferometer to L-band wavelengths (Howard Zebker, principal investigator) with the ad- 
dition of a second L-band antenna complementing the AIRSAR polarimetry aperture. The radar 
instrumentation was modified to allow operation in the usual TOPSAR mode, in which a single an- 
tenna transmits the radar pulses and both antennas receive, or in ping-pong mode, where alternate 
antennas transmit on alternate pulses but both antennas receive on all pulses. In ping-pong mode, 
the effective separation between the antennas, called the baseline is increased by a factor of two 
relative to TOPSAR mode, giving increased sensitivity to topographic height. This is particularly 
important at L-band, where the short baseline relative to the wavelength causes considerable height 
noise despite good signal-to-noise conditions. Typical imaging geometries and height sensitivities 
for TOPSAR are summarized in Table 1. 

The dual-wavelength, dual-baseline capability of the new TOPSAR system allows detailed study 
of the effects of vegetation canopies on the interferometrically derived heights. Coherent backscat- 
ter of microwave radiation is dominated by scatterers sized of the order of the radar wavelength. 
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It is well known from comparison of heights at field/forest stand boundaries in flat terrain that the 
C-band signal scattering in the forest is dominated by the canopy such that the interferometrically 
derived height is well above the forest base. The dual wavelength capability of the interferometer 
measures the differential penetration depth (if any) to provide insight into the scattering mech- 
anism, which can be complicated in real canopies. The vegetation canopy increases the effective 
scattering area of each resolution cell, leading to increasingly decorrelated radar echoes as the inter- 
ferometer baseline is increased. The dual baseline capability of the TOPSAR interferometer tests 
this prediction and establishes the usefulness of decorrelation measurements for estimating canopy 
depths. 

This work takes the first steps toward addressing these fundamental questions by first cali- 
brating the system at C- and L-band, then comparing C- and L-band height maps. Before the 
measurements can be interpreted, calibration is essential. The locations of the phase centers of the 
interferometer antennas are very difficult to measure at the millimeter scale required for accurate 
topographic mapping (Madsen et al 1993). Centimeter scale errors in baseline components can 
lead to planimetric and vertical map distortions in the tens of meters. In addition to the baseline, 
systematic phase errors exist in the radar data that are attributable to scattering of energy into an 
antenna from multiple paths, ostensibly the antenna tarings. These systematic phase errors cor- 
rupt the derived height maps in the cross track direction with vertical and horizontal oscillations 
of several meters’ amplitude at C-band and over 10 meters at L- band, severely limiting the utility 
of the measurements for sensitive differential penetration studies. 

Calibration 

The C-band TOPSAR baseline is measured yearly to an accuracy of a few centimeters. Until 
recently inaccuracies in the aircraft attitude and position measurement system precluded meter 
scale accuracy mapping. With the operational use of a refined inertial navigation system in late 
1994 and in 1995, a finer accuracy baseline for C-band is needed. Calibration using an array 
of corner reflectors in Rosamond dry lake bed is under way, however the flight conditions were 
turbulent and windy, so verifying the derived baseline parameters will require additional work. 
For the purposes of this study an alternate path of calibration was taken. It has the advantage 
of removing height distortions due to baseline errors and multipath phase errors simultaneously. 
This phase screen approach was first employed when JPL calibrated the Environmental Research 
Institute of Michigan’s IFSARE X-band topographic mapper (Kim and Imel 1995). 

The phase screen approach requires a reference digital elevation model. The TOPSAR digital 
terrain map is computed as an array of height values on a locally spherical grid. At the same time, 
the imaging geometry of each array element is recorded. The difference between the TOPSAR 
map and the reference DEM is taken - this constitutes an array of height errors. The residual 
interferometric phase error needed for calibration A<$ is related to these height errors Ah through 
the ambiguity height equation, given in a simplified two dimensional geometry by 

A , . . 27rn B cos (9 — a) . , . . 

(r = -r 0.1 

A r sin(0) 

where B is the baseline length, a is the angle the baseline makes with the horizontal measured 
counterclockwise, 6 is the look angle, A is the wavelength, and n is 1 for TOPSAR mode or 2 for 
ping-pong mode. The height-induced errors from both multipath and baseline error are primarily 


78 



functions of the look angle, so the phase screen correction inferred from the difference DEM and 
the geometry is tabulated as a function of look angle. By using the entire two dimensional DEM, 
the one dimensional phase screen benefits from averaging and is quite smooth despite errors in the 
DEM data. Using this phase screen to correct the interferogram, a new TOPSAR digital terrain 
map is produced and compared to the reference DEM. An adjusted phase screen is then computed 
and compared to the first estimate. This procedure is iterated until the phase screen becomes 
stable, usually after no more than two iterations. 

Figure 1 illustrates the phase distortions in the radar heights in the left-most panel. The area 
is located at Elkhorn Slough, CA and is quite flat, containing a large part of Monterey Bay in the 
left of the image. Despite the flat terrain, there is a tilt of the heights from left to right (near range 
to far range) in the image as well as oscillatory banding aligned with the flight direction, the effects 
of baseline error and multipath. The center panel shows the difference between the radar heights 
and the USGS 90m DEM resampled to the radar coordinates at 30 m spacing. The USGS DEM 
is the best readily available digital data set for this site. Nearly all the terrain features remaining 
in the difference map are due to inaccuracies or resolution limits in the USGS data. Nevertheless, 
the net phase screen derived from the difference map is quite smooth. 

On the second iteration, the region of the difference map used was restricted to the ocean and 
featureless areas in the bottom half of the image. This improved the estimate considerably. The 
right most panel of Fig. 1 shows the calibrated height data. The ripples in the Bay and flat regions 
clearly evident in the left panel have been removed in the right panel. If we assume that the 
average height of the USGS map over wide areas is correct, we can conclude that the phase screen 
calibration performed well. Figure 2 shows the phase screen profile as a function of look angle. 

The L-band system is new, so the baseline is commensurately less well known than the C-band 
baseline. To calibrate L-band, we decided to estimate its baseline by comparing the interferometric 
phase at C-band to the phase at L- band. If both baselines are perfectly well known then the 
following relationship holds: 


$ L (r) 


Ac B l cos (9(r) - a L ) 
X L B c cos (0(r) - a c ) 


$ c (r) = 0 


(0.2) 


Any residual phase ri<h(r) remaining as a function of range r we attribute at this stage to baseline 
error through the relation 

— c£$(r) = dB sin(0 — a) — daB cos (9 — a) (0.3) 

27rn 

It is straightforward to derive the corrections dB and da from the residuals d<&. The L-band baseline 
given in Table 2 was determined in this way. As a cross check of the baseline estimate, we measured 
planimetric distortion of the data before and after baseline correction by correlating the L- band 
interferogram with the C-band data. Before baseline correction, the images were misregistered by 
many pixels particularly in the near range. After baseline correction, the images aligned to within 
a pixel everywhere. 

The phase screen correction for L-band is still necessary because of multipath. Instead of using 
the USGS DEM, we used the newly calibrated C-band DEM. Few large scale artifacts exist in the 
difference map, however the large amount of height noise in the L-band data due to large ambiguity 
height limited the effectiveness of the phase screen. Residual multipath ripples remain in the L- 
band DEM at the 1-2 meter amplitude level. Though this is considerably lower than the statistical 
height noise at L-band, the artifacts are visible to the eye. 
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Preliminary penetration results 

Two differential penetration maps have been produced. The first is located at Elkhorn Slough, 
an obvious choice for calibration because the surface is essentially unvegetated. Figure 3 shows 
the calibrated C-band (left panel), L- band (center) and the C-L height difference (right panel). 
We expect the right panel to have no topographic variation even in areas of high relief, and this is 
observed: The L-band calibration succeeded, and the statistical height error of 2 - 5 m is sufficiently 
low to allow detection of penetration differences. 

Figure 4 shows the same three panels for a different site, Laurel Quad, CA located between San 
Jose and Santa Cruz. Laurel Quad is heavily vegetated with tall conifers and deciduous trees, as 
well as populated with grasslands and clear-cuts. Penetration differences should be apparent in the 
C-L difference map as positive numbers because penetration is expected to be larger at L- band. 
Indeed, close examination of the right panel shows virtually no height differences in the urban areas 
to the north and south but statistically significant differences in the forested areas between. Note 
that the imaging geometry and acquisition time for the C- and L-band systems are identical so 
there can be no terrain bias effect or angular calibration error causing this difference. Furthermore, 
careful study of the figure, or inspection of the full resolution data (these data are averaged for 
display), shows a clear correlation between backscatter variations associated with canopy and the 
height differences. In particular, above the image center southeast of San Jose, the northern slopes 
show clear vegetation, while the southern slopes axe bare and show no vegetation. 

Transects through the urban areas and forest consisting of profiles averaged across many image 
lines show clear distinctions, with an average difference of zero in the city and 4-6 meters in the 
forest. It is not clear how representative the average transect height is of the actual penetration 
locally. The image shows patchiness that is often attributable to variations in canopy. Because of 
the large height error in the L band data however, considerable averaging is required. Therefore, 
to better characterize the penetration, uniform tree stands must be identified from ortho images 
and other ground truth and a statistical average computed, probably as a single number per stand. 
This detailed analysis is under way. Comparison of these data with other acquisition modes will 
also provide insights into penetration phenomenology. However, calibration of the remaining modes 
is just beginning. 
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Table 1 

Relevant Interferometric Mapping Modes and Parameters^ 


Map Product posting (m) 


5 

Map Product Resolution (m) 


15 

Range Bandwidth (MHz) 


40 

Range Resolution (m) 


3.331 

Pulse Repetition Ground Spacing (m) 


0.34 

Number of Looks 


30 

Aircraft Altitude (km) 


8 

Range (km) 

9 

17.5 

Look Angle (deg) 

25 

63 

C/L-Band Ambiguity Height (m) 

56/243 

174/755 

C/L-Band Nominal Interferometric Correlation 

0.90/0.94 

0.90/0.94 

C/L-Band Nominal Statistical Height Error (m) 

0.56/1.81 

1.73/5.63 


* For 40 MHz ping-pong mode 


Table 2 

Preliminary TOPSAR Interferometric Calibration Parameters 


Wavelength 

L-Band (24 cm) 

C-band (5.6 

cm) 

Antenna Location 

Top 

Bottom 

T-B“ 

Top 

Bottom 

T-B 

X'f (m) 

0.54 

0.54 

0 

2.46 

2.46 

0 

Y(m) 

1.231 

0.52 

0.711 

1.154 

0.157 

0.997 

Z (m) 

1.601 

- 0.23 

1.831 

1.577 

-0.717 

2.294 

A <t>d (rad) ^ 

0.85 

0 

0.85 

0 

-1.35 

1.35 

A t d (ns) 

428.02 

444.72 

-16.7 

406.08 

437.04 

-30.96 


° T-B implies baseline components 


t X coordinate measured forward along aircraft fuselage 
Y coordinate measured outward to port side 
Z coordinate measured vertically upward 
Channel phase A <j>d and time delay A t c are mode dependent. 
These values are valid for 40 Mhz ping pong mode only. 
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TOPSAR "PHASE SCREEN" CALIBRATION 




Figure 1 


Elkhorn Slough, CA Run 180-1 1995 
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TOPSAR DUAL FREQUENCY TOPOGRAPHIC MAPPING 
Elkhorn Slough, CA Run 180-1 1995 
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L-BAND (24 cm X ) 
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Figure 3 
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TOPSAR DUAL FREQUENCY TOPOGRAPHIC MAPPING 
Laurel Quad, CA, Run 180-3 1995 
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MEASUREMENTS OF OCEAN FEATURES WITH AIRSAR 
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ABSTRACT 

Oceanographers have tried to obtain ocean data from SAR since mid 1978 
when SEASAT was launched and exposed an impressive array of features on 
the ocean surface. The difficulty in obtaining quantitative ocean data from SAR 
imagery is due to the fact that different ocean/atmosphere processes modulate 
the grey shades in the imagery. It is difficult to tell whether a front observed in 
the imagery is due to changes in the air/sea temperature difference, wind 
variability, horizontal current shear, rain front or spawning of near surface 
organisms. Determination of quantitative values, for the shear for example, are 
not realistic. With this gloomy introduction in mind, how do we proceed? 

Here we will discuss an indirect approach for measuring ocean features from 
SAR imagery. The approach is based on the refraction of swell that can be seen 
in the SAR imagery. Locally, swell refraction is due only to changes in 
bathymetry and surface current shear and it is easy to separate these two 
effects. Utilizing the indirect approach requires the following steps: 

1 . Image segmentation to a size smaller then the velocity scale. 

2. Dominant waves enhancement in each segment. 

3. Determination of the wave number of the dominant waves in each segment. 

4. An inverse model utilizing the refraction of the dominant waves. 

It has been shown already that linear energetic shears can be detected and 
measured from high resolution images of the ocean (see for example Sheres et 
al JGR, May 1985). A simple approach to measuring surface velocity and shear 
is to divide the imagery to smaller segments where the shear in each segment is 
approximately linear (not curved); the known linear algorithm would then be 
applied to each segment; the velocity /shear distribution can then be determined 
be assembling the segments with linear shear into an estimate of the actual 
shear/velocity distribution. 
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Swell is many times buried in the noise of SAR imagery. To use the swell in the 
SAR imagery as a diagnostic tool to measure shear/currents, we frequently need 
to enhance the dominant waves in the SAR imagery. An example of dominant 
wave enhancement using AIRSAR imagery will be presented. The imagery was 
obtained in May 26,1993 over a coastal area in the Gulf of Mexico exhibiting a 
rich array of features. 
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ABSTRACT 

The effects of terrain geometry result a great problem in SAR data calibration and interpretation 
of the relationships between geophysical parameters and radar backscattering measurements. The 
terrain correction could be done if we can determine the sensor position accurately and if high 
resolution and accurate Digital Elevation Model (DEM) is available at study site. It is clear that 
the radar backscattering highly depends on the local incidence angle. In this study, we have 
developed a model for the bare surface based on the simulations of Integral Equation Model 
(IEM) for a wide range of surface dielectric and roughness conditions. This model characterizes 
the relationship between the backscattering coefficients of W and HH polarizations at L-band 
and the incidence angle. The model does not require the surface dielectric and roughness 
information so that it is possible to determine local incidence angle for each pixel with only SAR 
measurements. We tested this model with NASA/JPL AIRSAR over Washita. The results showed 
an accuracy within a few degree can be achieved over the bare surface. However, the vegetated 
surface has significant effect on this relationship. 

An example of applications of this model is to estimate snow density. The imaginary part of 
dielectric constant of ice is very small so that the microwave is transparent over snow cover. The 
penetration depth could reach several ten meters. Since snow is a dense media, with density 
ranges from 0.1 to 0.55 in general, the real part of dielectric constant of snow is dependent only 
on the snow density. Our evaluation has shown the backscattering measurements at L-band over 
dry snow covered terrain are mainly affected by: 

1. Incidence angle change: from incidence angle at air/snow interface to refractive angle at 
snow/ground interface. 

2. The wavenumber (or wave propagation constant) will be shifted from k () to k 0 ^. 

£ £ 

3. The dielectric contrast will be changed from to ~ s ~. 

4. The transmission loss at the air/snow interface. 

Where subscript g, s, and 0 represent ground, snow, and air, respectively. With consideration of 
above factors, we can determine the refractive angle at sown/ground interface by using the model 
that we have developed. Furthermore, the snow density estimation can be obtained for each pixel 
by comparing the local incidence angle (at air/snow interface) calculated using DEM and the 
estimated refractive angle from L-band SAR measurements. We will demonstrate this technique 
and validations in detail. 
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1.0 Introduction 

Sandy barrier islands dominate the Atlantic and Gulf of Mexico shorelines of the United 
States. These islands naturally erode and accrete in response to waves, wind, tides, and long-term 
sea-level change. Since the late 1800s, barrier island shorelines on the whole have eroded at an 
average rate of 0.8 m/yr, although in many areas the rate can be several meters per year (Gibeaut 
and Crawford, 1994). The barrier island and peninsula complex off the coast of Texas is 
composed of several interactive depositional subenvironments. Tidal inlets (semipermanent 
breaches in the island) commonly open, close, and migrate along the coast with significant changes 
occurring over time scales of several months to decades. Storms can lower the beach and dunes by 
as much as 3 m while simultaneously depositing as much as a meter of washover sand in the 
backbarrier areas. They can also move shorelines tens of meters landward and create channels, 
scarps, and dune blowouts with a topographic relief of several meters. Post-storm recovery of the 
dune system, however, can take several years. During fair-weather conditions, the elevation and 
horizontal position of the beach changes daily from 1 cm to 1 m, depending on wave and wind 
conditions. Currently it is impossible to rapidly and accurately measure these changes across entire 
barrier islands or in various subenvironments. 

Galveston Island and Bolivar Peninsula reach a maximum elevation of 3 - 5 m above sea 
level (White et al., 1985). Because of their extremely low relief, even slight changes in the surface 
topography can result in significant variations in the soil and vegetation. Detailed knowledge of the 
topography can improve predictions for sediment erosion/accretion processes, storm impact and 
recovery, and the capability to map coastal wetland environments. Coastal geological studies have 
primarily focused either on long-term studies of large-scale phenomena at coarse spatial resolutions 
or on high-resolution, short-term studies of small-scale phenomena. Long-term studies are useful 
in describing general coastal change caused by cumulative effects of coastal processes. Short-term 
studies measure fine-scale changes in beaches and dunes occurring over periods ranging from one 
wave cycle to one monthly tidal cycle or during a single storm. They do not accurately show the 
integrated response of coastal depositional systems across a large area. To understand the 
dynamics of barrier islands, the interactions between various subenvironments and adjacent islands 
must be considered. A means to rapidly collect accurate high-resolution topographic data across 
tens of kilometers of shoreline has not existed, however. NASA/JPL's fully-polarimetric multi- 
frequency AIRS AR system has the potential to collect these data (using TOPSAR mode or repeat- 
pass interferometry techniques). When these topographic data are combined with wave, current, 
wind, and sediment data, a better understanding of the overall response of the coast to various 
environmental conditions will be possible. This understanding will in turn allow formulation of 
better predictive models of shoreline erosion. 

A three-year project to study small-scale topographic changes and relict geomorphic 
features on barrier islands using synthetic aperture radar (SAR) is underway at the University of 
Texas at Austin as part of a NASA program to study global topographic change. A study area on 
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the Texas coast consisting of Galveston Island and Bolivar Peninsula was overflown by the 
NASA/JPL DC-8 AIRSAR in April, 1995. Data was acquired in the fully polarimetric mode using 
C-, L-, and P-bands and in the TOPSAR configuration with C- and L-bands in interferometric 
mode. The study area will be overflown again in late spring of 1996. The data will be registered 
to GPS-surveyed points to form high-resolution digital elevation models (DEMs), and then 
analyzed to investigate possible topographic changes. 

Initially, the fully polarimetric data from a test site on Bolivar Peninsula is being used to 
model the SAR backscatter coefficient (a°) response to wetland vegetation. Marshes that 
commonly form the landward shorelines of barrier islands generally accrete and erode more 
slowly, but are often subject to rapid natural or man-induced changes. In the Galveston Bay area, 
about 20% of coastal marshes were lost from the 1950s to 1989. Draining, filling, and man- 
induced submergence caused by fluid withdrawal (water and oil) were responsible for most of the 
loss (Gibeaut and Crawford, 1994). By modeling the backscatter from the vegetation in the 
various sub environments, the study will investigate the use of SAR for studying wetland 
environments (Durden et al., 1995), (Ott et al., 1990). The vegetation species in this low-relief 
area are very sensitive to slight changes in the surface elevation. Using both GPS measurements 
and the TOPSAR data, the relationship between the vegetation and the topography will be 
investigated. 


2.0 Background 

Mathematical models have been developed that relate SAR backscatter measurements to 
surface parameters such as soil moisture (Oh et al., 1992), (Shi et al., 1992), (Dubois et al., 1995) 
and to vegetation parameters like woody biomass and branching structure (Dobson et al., 1992), 
(Durden et al., 1989). However, most of the previous research in these two areas has been 
devoted to the study of soil properties for bare to sparsely- vegetated terrain and vegetation 
properties for agricultural crops and forested areas. Studies have also demonstrated the capability 
to distinguish between wetland sub environments using SAR, (Hess et al., 1995) and (Pope et al., 

1994) , and more recently, the ability of vegetation models to explain the observed backscatter from 
herbaceous vegetation in the presence of surface water (Durden et al., 1995). Many of the models 
developed to represent the o° response to natural terrain either describe scattering from a surface or 
scattering from vegetation above the surface. Each class of model relates the observed a° to 
parameters that characterize the geometric and electrical properties of the study area. Once a 
mathematical relationship has been developed, parameters of interest are computed via model 
inversion. 

Surface models such as (Shi et al., 1992) and (Dubois et al., 1995) attempt to represent the 
different polarization combinations of o°, and their linear combinations, as functions of the angle of 
incidence, geometric properties ( e.g. rms surface height), and electrical properties ( e.g. relative 
dielectric constant ). Because of the large difference in electrical conductivity between water and 
diy soil, the relative dielectric constant can be used to estimate the soil moisture using dielectric 
mixing models (Dobson et al., 1985). The mathematical expressions relating a 0 to these 
parameters are either based upon electromagnetic theory or empirically derived (Chen and Fung, 

1995) . Theoretical models such as the small perturbation method and the phase perturbation 
method are useful for investigating the dependence of a° on various parameters, but are only valid 
for surfaces that are very smooth relative to the incident wavelength (Oh et al., 1992). Most 
natural surfaces exceed this limit on surface roughness. Although the integral equation method, 
which may apply to a greater variety of surfaces, has been recently developed (Fung et al., 1992), 
researchers have historically utilized empirical surface models to explain o° behavior over natural 
terrain. 
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Empirically-based models can better 'deal with the levels of surface roughness encountered 
in natural terrain, but they suffer from two major shortcomings: (1) a strong dependence upon the 
data used to develop them, and (2) difficulty in separating the effects of surface roughness from 
those of soil moisture. Once a backscatter model has been developed by empirically fitting a 
function to a particular data set, it usually does a poor job of predicting o° for other data sets 
(Dubois et al., 1995). Although this problem can be alleviated by using several data sets to 
develop the model, it will still perform poorly if the terrain being analyzed differs significantly 
from that used to develop the model. The second shortcoming of the empirical surface models is 
their difficulty in separating the contributions of surface roughness and dielectric constant relative 
to the total c° in the presence of significant vegetation. The primary effect of the vegetation layer 
covering the surface is to introduce a volume-scattering component into o°. This produces the 
same response that an increase in either surface roughness or the relative dielectric constant 
(decrease in soil moisture) would cause, so the change in o° cannot be uniquely attributed to either 
effect. As a result, these models tend to overestimate surface roughness and underestimate soil 
moisture when vegetation is present. 

The impact of vegetation cover on the predictive capability of a model varies substantially. 
Researchers have investigated various approaches to reduce the sensitivity of the models to 
vegetation cover. Ratios of the co-polarized o° responses and their linear combinations have been 
considered (Dubois et al., 1995) as well as accounting for the volume scattering and surface 
scattering components of o° separately (Shi et al., 1992). These approaches often achieve 
improved representation of the measurements, but the problem has not been completely solved, so 
models should only be applied to appropriate data. The Dubois et al. (1995) model includes a 
criterion for determining whether the vegetation cover is too thick for the model to provide 
adequate results. Areas where this criterion is exceeded are excluded from consideration. The 
a°h v /a°vv ratio is used as an indicator of significant vegetation since the amount of diffuse (multiple) 
scattering usually increases in the presence of vegetation. A threshold of -1 1 dB at L-band was 
determined empirically for this particular model. Applying this criterion to the Bolivar Peninsula 
test site excluded approximately 85% of the terrain. Because this class of models is not directly 
applicable to the data acquired of the Bolivar test site, and appropriate generalizations are not 
apparent, models of the scattering due to vegetation are being investigated instead. These models 
will be used for discriminating vegetation in the test site and representing the behavior of o° in the 
presence of marsh vegetation. 

Vegetation models have been used to represent the plant material both as a layer of 
continuous and discrete random media (Lang, 1981). Several techniques for analyzing the 
scattering from continuous random media were developed in the 1960s and 1970s. These models 
assume that the permittivity of the medium is a random variable whose statistical moments are 
known, and then determine the average backscatter cross section from the layer. More recently, 
the plant material has been treated as a layer of randomly-oriented discrete scattering elements. 

This approach allows a more direct representation of various plant elements ( e.g. trunks and size, 
shape, and orientation of leaves and branches). These plant elements are often described by 
dielectric cylinders and disks with random orientations prescribed by probability distribution 
functions (Yueh et al., 1992). Much of the research with discrete vegetation models has focused 
on agricultural crops and forested areas (Brakke et al., 1981) and (Durden et al., 1989). Recently, 
these types of models have also been applied to grasslands and wetlands (Saatchi et al., 1994) and 
(Durden et al., 1995). In particular, Durden et al. (1995) focused on flooded herbaceous 
vegetation similar to some areas of the Bolivar Peninsula test site. 

These vegetation models can also be considered to be empirically-based in that model 
coefficients and characteristics of the scattering elements can be varied to achieve a better fit to the 
observed o°. Insight can also be gained into the scattering behavior of the vegetation by varying 
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the orientation, size, and shape of the scattering elements. Some of these models can also be used 
to determine what fractions of the total c° are due to various scattering mechanisms, including 
direct scattering from the vegetation and direct-reflected scattering from ground/vegetation 
interaction (Lang and Sidhu, 1983). The model from (Durden et al., 1995) expressed the 
attenuation coefficient and scattering cross-section per unit volume in terms of the leaf area index 
(LAI). It was then possible to estimate LAI, which can be a useful parameter for assessing plant 
growth stage and density. 


3.0 Site Description and Data Collection 

The initial test site selected for developing the wetland-vegetation model is a saltwater 
marsh on the landward side of Bolivar Peninsula, TX (Fig. 1). A single fully polarimetric, three 
frequency frame from the AIRSAR coverage over the Bolivar Island test site is being analyzed 
initially (Fig. 2). The site contains a variety of terrain types and associated vegetation. The upland 
flats of the peninsula consist mostly of herbaceous range grasses. Towards the mainland, the 
terrain gradually decreases in elevation and supports a saltmarsh. Between the uplands and 
saltmarsh lies a highly saline intertidal zone where salt flats and salt barrens form. The saltmarsh 
can be further segmented into lower and upper marshes which correspond to frequently flooded 
and non-flooded areas. The elevation ranges from sea level at the edge of the saltmarsh to about 2 
m on the upland flats, with the highest point occurring along the foredunes at approximately 4 m. 

The vegetation in this test site is strongly dependent upon the elevation. The saltmarsh 
becomes less frequently flooded further from the landward side of the peninsula and the tidal 
inlets. The intertidal zone develops at the maximum elevation where flooding and soil moisture are 
still due to seawater. Here, the salt compounds are precipitated out of the water through 
evaporation, creating a highly saline environment. Finally, on the upland flats, the elevation is 
sufficient to prevent the frequent flooding. Each of the regions can be distinguished by the most 
commonly occurring vegetation and, to a lesser extent, the soil type. For developing the vegetation 
models, the test site is considered to be comprised of upland flats, an intertidal zone, non-flooded 
saltmarsh, and flooded saltmarsh. The upland flats, near the center of the peninsula, are 
characterized by loose, sandy soil. The water table is often close to the surface, but the topsoil is 
still typically well-drained and dry. This region is covered with range grasses such as Evagrostis 
(Webb), which tend to be short (8-10 cm) and have blades oriented between 45° and 90°. The salt 
barrens and flats of the intertidal zone also have sandy soil, but are more saline. This region 
supports only sparse vegetation such as the short grasses monathochloe littoralis, and succulents 
such as Salicornia (Webb) and (White and Paine, 1992). The monathochloe littoralis has very 
short leaves (5-10 mm) which occur along the ground runners. The Salicornia grows erect to 
heights of 6 cm. The high, less-frequently flooded saltmarsh occurs just below the intertidal zone. 
Marsh soils have a higher concentration of silt and organic material than the upland terrain. The 
high marsh is dominated by tall grasses called spartina patens, which grow to 90 cm and are 
vertically oriented. The lower saltmarshes are very often submerged under several centimeters of 
standing water. The dominant species is a very tall grass (150 cm) known as spartina altiniflora. 

As expected, preliminary visual interpretation of the image showed that C-band is the best 
frequency for differentiating vegetation differences. The longer wavelength L- and P-bands 
penetrated the vegetation and revealed large scale surface variation which is in rough agreement 
with the soil variation on the peninsula. The sandy soils of current and past beach ridges produced 
stronger a° responses than the silty soils that comprise the existing marsh over which the beach 
ridges have accreted. Field measurements will be used to determine whether the radar is actually 
penetrating the vegetation and sensing soil differences or merely indicative of variation in the 
vegetation that itself is responding to soil differences. 
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Combinations of various polarization's were employed as simple indicators of the relative 
abundance of vegetation. The difference between the P-band co-polarized data (c 0 ^ - a° vv ) was 
examined because it was expected that the vegetation would attenuate a° vv more than o 0 ^, 
especially at longer wavelengths. The response was relatively weak throughout the intertidal zone, 
which agreed with the sparse vegetation in that region. Ratios of cross-polarized values of 0 0 were 
also investigated. The L-band ratio 0 °h v /c° vv was expected to indicate the presence of 
depolarization caused by multiple scattering from very rough surfaces and vegetation, with a 
reduced dependence upon local slope and incidence angle. The response was weak in the intertidal 
zone and very strong in the uplands, with an intermediate level of response over the marsh. The 
intertidal zone produced the weakest response because of the sparse vegetation. The relatively 
weak response over the marsh is believed to be due to the specular nature of the standing water and 
extremely smooth surface in that region. 


4.0 Analysis & Future Work 

Ground measurements will be acquired from these four regions to support the development 
of the vegetation models. Surface parameters such as rms surface roughness, volumetric soil 
moisture, and saline concentration will be measured. Dielectric mixing models will be used to then 
estimate soil dielectric constant. (Dobson et al., 1985). Plant geometry and moisture content will 
also be measured. Average plant height, size, shape of the leaves, and leaf and stalk orientation 
will represent the size and orientation of the scattering elements in the models (Lang, 1981). Plant 
moisture content will be used to estimate vegetation dielectric constants (El-Rayes and Ulaby, 
1987). Adjustment will be made for differences in plant geometry due to annual growth which has 
occurred since the SAR data was acquired. Model results will be compared to the o° data from an 
AIRS AR frame centered over the test site. Adequate geometric correction of the AIRS AR frame 
will be achieved using extensive GPS data and GPS-surveyed comer reflectors visible in the 
image. 


The structural form of the model developed by (Lang and Sidhu, 1983) is being used 
initially. Vegetation is treated as a layer containing scattering elements over a flat lossy ground, or 
a specular reflector if flooded (Durden et al., 1995). The herbaceous vegetation, which mostly 
consists of grass stalks and blades of varying height and orientation, is modeled as a layer of 
dielectric disks with specified orientation distributions. The dielectric constant of the ground is 
being adjusted to account for the moisture and salinity variations observed in the field 
measurements. The model uses the Foldy technique to determine an approximate equation for the 
mean field. Using this mean field, an equivalent dielectric constant for the scattering layer is 
determined. The correlation of the field is then determined by assuming that the scatterers are 
imbedded in the equivalent medium (i.e. the distorted Bom approximation). The model results are 
being compared to the SAR data to investigate differences in the predictions of the a° values for the 
various subenvironments. 

The configuration of the scattering elements will be investigated to best approximate the 
true plant geometry. Other wetland sites in the study area will also be analyzed to test the range of 
applicability of the model. Finally, the correspondence between vegetation type and topography 
will be examined. 
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Figure 2 Test site ( C-band total-power image acquired by the NASA/JPL AIRSAR) 
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Enhancements of AIRSAR data have demonstrated the benefits of radar for revealing an 
additional and much higher level of information about the composition of the terrain than 
enhancements of either SPOT-PAN or Landsat TM data. With the help of appropriate 
image-processing technology, surface and near-surface geological structures, hydrological 
systems (both current and ancient) and landform features, critical to the comprehension of 
adjacent mineralisation but often obscured by vegetation and unconsolidated materials, have 
been put plainly in evidence in a diverse range of landscapes. 

In the Great Sandy Desert region where spectral variability is minimal (as observed by 
Landsat TM), radar’s sensitivity to the micromorphology of sparse exposures of subcrop 
and lag gravels has provided a new insight into the regioms geological framework, its 
landforms and their evolution that is unobtainable through conventional methods. Subtle 
variations in surface roughness of terrain components at scales near the observing 
wavelengths is the dominant factor contributing to the radar backscatter from the landscape 
components. Since most of the GSD is masked by a layer of sand and vegetation, and that 
much of the rock exposure is intermittent, localised and scattered, the unique ability of 
enhancements of AIRSAR data to highlight new details of this arid landscape is striking 
when compared to the published geological map and the SPOT-PAN and Landsat TM 
scenes of the region. 

In the Pilbara region, advanced processing of AIRSAR data to unmix the backscatter 
between and within the three frequencies of data has highlighted sub-surface extensions of 
greenstone lithologies below sand cover and morphological evidence of past flow 
conditions under former climate regimes. The identification of palaeolandforms, subtle 
hidden geological structures, sequences of valley-fill sediments, regolith-landforms and 
sequences of landform development have all been achieved from enhancements of AIRSAR 
data. 

On the basis of these observations which are described fully in restricted CSIRO- AMIRA 
reports, it is recommended that radar remote sensing technology involving the use of high 
resolution, polarimetric data be seriously considered as a viable tool for exploration in 
erosional and depositional environments located within Australials mineral and oil- 
prospective provinces. 
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In an exploration program it is important to understand the regolith-landform relationships so as to 
allow appropriate geochemical sampling methods to be chosen, to permit sampling programs to be 
optimised and to permit sensible data presentation and interpretation. This presentation discusses 
the distinctive contribution AIRSAR data can make in characterising the regolith-landforms in a 
relatively vegetation-free environment 

AIRSAR frame-processed data were initially MAF-cleaned to enhance the signal content of the data 
before geocoding to AMG coordinates. Colours in a three frequency, single polarization 
combination image CW/LW/PW, and in an enhancement of pedestal height, relate directly to the 
scale of surface roughness of the various regolith units. Consequently, the AIRSAR images were 
able to provide excellent discrimination between and within regolith-landform units mapped in a 
previous study using air photos and Landsat TM data. 

Examination of the AIRSAR enhancements reveals that in mafic terrains, and to a lesser extent, in 
felsic terrains, AIRSAR data provides discrimination between the principal geomorphic regimes, 
relict, erosional and depositional. An important observation is the ability to consistently interpret 
the level of stripping of the landscape within the erosional regime from which details of regolith 
stratigraphy can be inferred. The multi-frequency signals respond differently to the relative 
abundance and texture of the surficial materials including soil, lag gravels and outcrop, for each 
unit within the weathering profile. Whilst the C-band response is similar (high) for all units, L- 
and P-band signals vary significantly from low to high according to the presence and texture of the 
surface float and outcrop. 

A multivariate statistical technique called an all-possible subsets calculation was used to examine 
the degree of polarimetric separation between selected regolith-landforms for all combinations of 
the nine-band AIRSAR data. This process uses, as a statistical base, the pixel values of training 
sites representing the regolith-landforms which are important for the polarimetric separation of the 
regolith-landform classes. In this study, the process indicated the most significant 3-band 
subset(s) suitable for enhancing the polarimetric contrast between regolith-landforms on an image 
processing screen or in a hardcopy image. 

An unanticipated aspect of the research has been has been the identification on the AIRSAR 
imagery of previously unmapped, structural features, these same features cannot be identified in 
air photos or other remotely sensed datasets. Their mapping has led to the development of a 
hypothesis for describing the structural control and development of landforms within this region. 

The requirement for a regolith-landform framework and knowledge of the processes of landscape 
evolution is integral to the planning, execution and interpreting of geochemical sampling programs. 
Of necessity this entails the recognition and characterisation of regolith units and regolith 
stratigraphy, and the ability to distinguish between relict, erosional and depositional materials. 
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Abstract — SAR polarimetric processing has received a great deal of attention for the purposes of 
enhancing the detection, classification and/ or identification of scene scattering features. This paper 
describes and compares the results obtained from applying the polarimetric processing techniques that 
have been developed and/or advanced by researchers at MIT Lincoln Laboratory, and basic polarimetric 
matched filter (PMF) techniques that have been modified to treat distributed scatterers. The paper 
describes the individual polarimetric processing techniques and the formulations used for their 
application to polarimetric imagery obtained from the NAWC P-3 SAR. For all cases the polarimetric 
techniques have been applied to enhance the detection of distributed scatterers in clutter. The motivation 
for the emphasis arises from considering most complex scatterers of interest (either man-made or natural) 
as distributed scatterers consisting of multiple scattering centers, and many modem polarimetric radar 
systems have the resolution performance to resolve the multiple scattering centers. Furthermore, the 
individual scattering centers of a distributed target can often exhibit different polarimetric scattering 
characteristics, and consequently do not respond favorably to polarimetric processing techniques that 
have been derived for individual point scatterers. The treatment of distributed scatterers also includes 
concepts for tuning the polarimetric ensemble response of individual scattering centers. 
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1. Introduction 

While interferometric synthetic aperture radar (INSAR) was primarily developed as a tech- 
nique for measuring surface topography [Graham, 1974], it is also very sensitive to the 
vertical structure of vegetation. Parameters describing the vertical distribution of vegeta- 
tion, the simplest of which is the depth of the vegetation, are of interest for vegetation 
monitoring and ecological successional modeling. For example, vegetation layer depth and 
scattering characteristics can be used to infer biomass and age of forest stands [Manning 
et al. 1984]. Unlike single-aperture synthetic aperture radar (SAR), INSAR is sensitive to 
the distribution of the altitudes of the returns from a vegetation layer. The altitude distri- 
bution of returns depends on both the actual physical vertical distribution of vegetation in 
the layer, as well as the attenuation level of the radar signal at different depths in the layer. 
INSAR is therefore primarily sensitive to the physical structure of vegetation layers; it is 
affected by the chemical structure only through the degree to which attenuation, mainly 
due to water content, weights the returns from different depths into the vegetation. Because 
of its dominant sensitivity to physical structure, INSAR is potentially complementary to 
infrared and optical remote sensing techniques which axe primarily sensitive to chemical 
composition [e.g. Wessman et al. 1988], but can also be used to infer structural changes 
in vegetation layers [e.g. Vogelmann et al. 1989]. Because radar penetrates the vegetation, 
INSAR is also sensitive to the surface topography underlying vegetation layers, which is 
important for hydrology and the study of lateral vegetation distributions. 

Drawing from recently submitted work [Treuhaft et al., 1995], this paper first gives a 
heuristic description of the sensitivity of INSAR to vertical vegetation distributions and 
underlying surface topography. A parameter estimation scenario is then described in which 
the INSAR cross-correlation amplitude and phase are the observations from which vege- 
tation and surface topographic parameters are estimated. It is shown that, even in the 
homogeneous-layer model of the vegetation, the number of parameters needed to describe 
the vegetation and underlying topography exceeds the number of INSAR observations for 
single-baseline, single-frequency, single-incidence-angle, single-polarization INSAR. Using 
ancillary ground-truth data to compensate for the underdetermination of the parameters, 
forest depths are estimated from the INSAR data. A recently-analyzed multibaseline data 
set is also discussed and the potential for stand-alone INSAR parameter estimation is as- 
sessed. The potential of combining the information content of INSAR with that of in- 
frared/optical remote sensing data is briefly discussed. 
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2. Modeling the INSAR Response to Vegetation Vertical Distributions and 
Underlying Surface Topography 

INSAR consists of transmitting a radio signal from one or both of two antennas and 

receiving the signal scattered from the surface and/or vegetation at both antennas. The 

— ♦ 

antennas are located at either end of a baseline B, as shown in Figure 1. The Figure shows a 
horizontal baseline and horizontal ground surface for simplicity, with radar incidence angle 
do. The two receivers and (at least) one transmitter can either be mounted on an aircraft, 
as in TOPSAR [Zebker et al., 1992], or on dual spacecraft, as in proposed orbiting interfer- 
ometric missions. The vegetation is pictured as starting at the ground surface altitude of 
zq and extending up to zo + h v , where h v will be called the vegetation layer depth. Returns 
from the lower edge of the vegetation at Ro and from a vegetation constituent at Rj are 
shown. 



The fundamental INSAR observations are the amplitude and phase of the complex cross- 
correlation of the signals from the 1 and 2 ends of the baseline. As shown in Figure 2, 
the qualitative effect of vegetation on the magnitude of the cross-correlation is to reduce it 
relative to that induced by bare surfaces. This is because each scatterer at Rj with altitude 
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Imaginary Part of the Cross-Correlation 



Figure 2: The schematic effect of vegetation at various altitudes on the real and 
imaginary parts of the INSAR cross-correlation function. 

Zj contributes a phasor to the cross correlation proportional to exp(i<f>(zj,rij)) where 


= ko(rij — r 2 j)| 0 ~ k 0 B sin [cos 1 C 1 ) 


where ko is the wave number at the center of the bandpass, h is the altitude of the radar, 
and r\j is the range from the 1 end of the baseline to the scatterer j and similarly for r%j ■ 
The incoherent addition of the phasors at different altitudes will cause the aggregate cross- 
correlation to be of smaller amplitude than the zero-baseline case, shown along the abscissa 
in Figure 2. The qualitative effect on the phase, due to the addition of cross-correlation 
phasors, is to increase it to some value between that observed for bare surfaces and that 
observed for a surface at the top of the vegetation layer. Expressions for the complex cross- 
correlation are derived in Treuhaft et ah, 1995, on which most of this discussion is based. 
The expression for the INSAR cross-correlation in the presence of a vegetation layer with 
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vertical structure is given in an appendix in that reference. It is: 
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where A is the 1/distance loss due to spherical wave propagation, ko is the wave number 
of the radiation at the center of the bandpass, (ri — r 2 )|o is the difference in path length, 
between the two ends of the baseline, for signals originating from the point Rq in Figure 
1. In (2), W v is an instrumental azimuthal resolution function, with azimuthal coordinate 
T] (in and out of the paper in Figure 1), W r is an instrumental range resolution function, 
ro is the distance from the transmitter 1 to the point Rq in Figure 1, a r is the derivative 
of interferometric phase with respect to the range direction, and a z is the derivative of 
interferometric phase with respect to the vertical z direction. Both a r and a z depend on 
the baseline length B , the incidence angle, and the frequency of radiation. The vegetation 
backward scattering amplitude /t, the forward scattering amplitude ff, and the vegetation 
density p are all described as functions of altitude. The integral inside the square brackets 
at the end of (2) describes the attenuation in the medium. The phase in the argument of 
the first exponential is that typically associated with bare-surface INSAR and is given by 


0o = 4>(z 0 , r 0 ) = k 0 (ri — r 2 )| 0 ~ k 0 B sin [cos 



(3) 


where zq is the altitude of the underlying ground surface. The ensemble average brackets 
in (2) < • > represent averages over different configurations of scatterers and different 
scatterer compositions and orientations. From (2), it can be inferred that although INSAR 
is sensitive to the vertical distributions in vegetation density and scattering amplitudes 
(which depend on both scatterer geometry and chemistry), the amplitude and phase of the 
cross-correlation < E(Ri)I?*(-R 2 ) > alone are not enough to invert (2) for < / 2 (z) >, p(z), 
and Im < ff(z) >. 


3. Parameter Estimation Scenarios and Demonstrations 


Ultimately, INSAR data from many baselines and frequencies could be used to invert (2) 
for detailed vertical structure information, but for the limited baselines and frequencies 
available with TOPSAR, a homogeneous vegetation model is adopted, in which the functions 
of z in (2) are assumed to be spatially constant. For the homogeneous model, the complex 
cross-correlation becomes: 
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where </ 2 >, Jm<//>, and p 0 are now taken to be constant. For current baseline- 
frequency combinations, the homogeneous model is expected to work well for many vegeta- 
tion layers, as indicated in the reference. Now the three unknown parameters in Equation 
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(4) are 1) h v the vertical extension or depth of the vegetation, 2) zq the altitude of the 
ground surface (taken to be the lower edge of the homogeneous vegetation layer), and 3) 
a x the extinction coefficient, defined as 

4irp 0 Im < ff > 

cr x = ; 

Kq 

and is a measure of the power loss per unit length into the vegetation layer. 

The amplitude of the normalized cross-correlation, which is 


( 5 ) 
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due to vegetation is shown in Figure 3, for the Seven Mile Island area of the Bonanza 
Creek Experimental Forest (BCEF). Effects due to noise and correlation loss over the range 
resolution cell have been removed, as described in the reference. The areas marked “CC” in 
Figure 3 are clear cuts, those marked “WS” are white spruce, and those marked “BP” are 
balsam poplar. From Figure 2 and the discussion following it, it is expected that the cross- 
correlation amplitude should decrease (be blacker in the Figure) with increasing vegetation 
height. It can be seen that the clear cuts are whiter than the vegetated areas. In order to 
assess the sensitivity of INSAR to vegetation depth differences, note that the stand WS23 
(ground-truth heights 12m) is clearly brighter than WS5 (ground-truth height;«21m). It is 
fair to expect that parameter estimation of h v may ultimately exhibit few-meter accuracy. 

In an analysis in Treuhaft et ah, 1995, which used ancillary extinction coefficient data to 
compensate for the underdetermination of parameters mentioned above, forest layer depths 
h v were estimated from the BCEF data. The large regions, a few hundred meters on a side, 
shown in Figure 3 were used to estimate vegetation layer depths. Figure 4 shows the h v re- 
sults versus ground-truth tree height for two sets of extinction coefficients, one appropriate 
for frozen conditions and one for thawed conditions. Although the INSAR data were taken 
in the summer of 1993, and h v estimates should have been more accurate using thawed a x 
values, the opposite is true. Several reasons for the improved results (closer to the y=x 
line) with the frozen a x 's are given in the reference. These reasons include instrumental 
correlation losses, which would cause h v estimates to be too large, and the tendency of 
vegetation gaps to lower the effective extinction coefficient. Nonetheless, Figure 4 suggests 
that few-meter vegetation depth estimation may be possible in the near future. Sensitivity 
analyses suggest that the same accuracy can be obtained for the underlying surface topogra- 
phy, and that 0.1 db/m accuracy may be achievable for the extinction coefficient estimates, 
for sufficiently tall vegetation. 

In another demonstration which is not yet fully analyzed, multibaseline TOPSAR data 
were acquired by using the so-called “ping-pong” mode. It can be shown that, while the 
usual TOPSAR 2.5-m baseline is obtained when a signal is transmitted at one end of the 
baseline and received at the 1 and 2 ends, an effective 5-m baseline can be obtained if trans- 
mission and reception occurs at both ends of the baseline. Multibaseline, multifrequency 
(L and C band) data from the BOREAS Southern Test Site near Prince Albert, Canada, 
and the Kellogg site in Michigan were acquired and processed. Preliminary results indicate 
that the three parameters h v , zq, and a x can be estimated from the multibaseline INSAR 
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Figure 3: The normalized cross-correlation for the 7-Mile Island Area of the 
BCEF, corrected for correlation loss due to noise and finite range resolution 


data alone, but only for very steep incidence angles< 30°, and there appear to be systematic 
errors at the 10-m level for h v currently under investigation. Longer baseline pairs, for 
example 5 and 10-m baselines, may be better suited to accurate parameter estimation. The 
terrain selected for both demonstrations was flat, in order to avoid the complication of 
removing the correlation effects induced by surface slopes [Rodriguez and Martin, 1992] 
before estimating vegetation and surface parameters. 
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4. Combining the Information Content of INSAR and Infrared/ Optical Data 

INSAR and infrared (IR)/optical data have complementary sensitivities to vegetation char- 
acteristics, and may therefore be profitably combined in vegetation parameter estimation. 
For example, for shallow ( h v < 10m) vegetation layers, the INSAR determination of the 
extinction coefficient and h v is poor. This is because INSAR relies on the vertical diversity 
of returns to determine correlation amplitude and phase, from which the two parameters are 
inferred. When the vertical diversity is restricted for shallow layers, the extinction coefficient 
and the h v parameters become highly correlated and their accuracy can be compromised. 
IR and optical techniques are also sensitive to the extinction coefficient [e.g. Myneni et 
al., 1989], and it is possible that the extinction coefficient could be inferred from spectral 
rather than geometric signatures. In that case, it may be possible to use extinction coef- 
ficient information from the IR/optical data to infer vegetation densities through (5), and 
also appropriately transform the imaginary part of the average scattering coefficient to ra- 
dio frequencies. If it were possible to determine the radio extinction coefficieint from the 
IR/optical data as suggested above, introduction of the radio extinction coefficient into the 
INSAR analysis would enable much more accurate determination of vegetation depth, h v . 
This is an example of a potential combination of INSAR and IR/optical data which would 
yield more than either data type could supply separately. 
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5. Summary 

Homogeneous vegetation layer modeling enables the description of the INSAR response to 
vegetation in terms of three parameters, the layer depth, the extinction coefficient, and the 
altitude of the underlying ground surface. With many baselines and frequencies, modeling 
the INSAR cross-correlation could enable more-detailed extraction of the vertical profiles of 
vegetation. With current TOPSAR data over the BCEF, a demonstration using ancillary 
extinction coefficient data showed few to 10 meter agreement in the estimated vegetation 
layer depth. Another demonstration using multibaseline, multifrequency TOPSAR data 
over the Boreal forests and over the Kellogg site is in progress. In the future, combining 
INSAR observations with infrared and optical remote sensing in a unified parameter esti- 
mation scenario may yield much more complete information about vegetation layers and 
underlying topography. 
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Abstract: "WISIP: Wideband (PHz - PHz) Interferometric Sensing and 

Imaging Polar imetry* has become an important, indispensable tool in 
wide area military surveillance and global enivronmental monitoring of 
the terrestrial and planetary covers. It enables dynamic, real-time 
optimal feature extraction of significant characteristics of desirable 
targets and/or target sections with simultaneous suppression of 
undesirable background clutter and propagation path speckle at 
hitherto unknown clarity and never before achieved quality. "WISIP' 
may be adopted to the Detection, Recognition and Identification (DEI) 
of any stationary, moving or vibrating target or distributed scatterer 
segments versus arbitrary stationary, dynamically changing and/or 
moving geo-physical/ecological enivronments , provided the 
instantaneous 2^2 phasor (Jones/Sinclair) and 4x4 power density 
(Mueller/Kennaugh) matrices for forward-propagation/backward- 
scattering, respectively, can be measured with sufficient accuracy. 

For example , the DEI of stealthy, dynamically moving and/or 
camouflaged stationary objects occluded deeply into heterogeneous 
stationary and/or dynamically moving inhomogeneous volumetric scatter 
environments such as precipitation scatter, the ocean sea/ lake surface 
boundary layers, the littoral coastal surf zones, pack-ice and snow or 
vegetative canopies, dry sands and soils, etc., can now be 
successfully realized. A comprehensive overview is presented on how 
these modern high resolution/precision, complete polarimetric 
coregistered signature sensing and imaging techniques, complemented by 
full integration of novel navigational electronic tools, such as DGPS, 
will advance electromagnetic vector wave sensing and imaging towards 
the limits of physical realizability. Various examples utilizing most 
recent image data take sets of the NAWC/IEEIM-P3-UWB-T0PIF 5 E-CATT/LTBL- 
POLSAR and NASA-JPL-AIRSAR airborne , the NASA/DARA/I»SI-SIR-C/X-SAR 
shuttle, and the ESA ERS-1/2 satellite imaging systems will be 
presented for demonstrating the utility of WISIP. 
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l. Background 


GeoSAR is a project to develop a commercial airborne, radar-based, terrain mapping system 
utilizing technology originally developed by NASA and ARPA for commercial and defense applications. 
The project team consists of the Sensor Technology Office within the Federal Government’s Advanced 
Research Projects Agency (ARPA), NASA’s Jet Propulsion Laboratory (JPL) at the California Institute 
of Technology, Calgis Inc. a California business, and the California Department of Conservation (DOC) 
which brought together the GeoSAR team. GeoSAR has three primary goals. 


1. Validation of radar terrain mapping technology developed by NASA- JPL and ARPA. 


2. Application of this technology to seismic hazards and environmental mapping for the 
Department of Conservation and other governmental agencies. 


3. 


Demonstration of commercial marketability via end-user projects using real world problems and 
applications. 

2. Multi-Year Project 


The work plan is a 4 year end-user driven application and systems development and verification 
process, with the final goal being commercialization of the technology to serve the end-user community. 
Year 1 efforts focused on accuracy assessment for both terrain modeling and land cover analysis. Year 2 
efforts are now focusing on development of foliage penetrating (FOPEN) capabilities for mapping true 
ground surface and expanding end-user applications. Years 3 and 4 will concentrate on systems’ 
development and commercialization. 


3. YearI- Overview of Activities 


Year 1 began in the spring of 1994 and focused on end-user experiments designed to test and 
validate interferometric synthetic aperture radar (IFSAR) mapping capabilities using data acquired 
from JPL’s airborne TOPSAR (Topographic Mapping Synthetic Aperture Radar) instrument flown 
aboard the NASA/JPL DC-8 research aircraft. 


A ground campaign was conducted to establish geodetic control, set radar reflectors, and acquire 
ground truth for four test sites. The sites were chosen to represent environmental mapping problems 
from DOC programs that could be solved using a GIS (geographic information system) enhanced with 
radar mapping technology. JPL oversaw flight planning, acquired, and processed the TOPSAR radar 
data. Calgis led the ground campaign and provided end-user support. DOC provided end-user input, 
and conducted application and verification tests. 


* This paper is in part excerpted from the GeoSAR Year 1 Research and Development Status Report. 
prepared for the Advanced Research Projects Agency, by the: California Department of Conservation, 
Sacramento, California; the Radar Sciences Section, Jet Propulsion Laboratory, Pasadena, California; 
and, Calgis Inc. Fresno, California. 
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Between May and July 1994, JPL conducted three TOPSAR flights that acquired radar data over 
five sites. 

• Antioch Hills in the Sacramento Delta area serves as a super site for three end-user 
experiments continuing into Year 2 : (1) digital elevation model (DEM) verification; (2) slope 
hazards assessment; and, (3) land use and land cover analysis. 

• Laurel Quad, north of the town of Santa Cruz; is the primary site for foliated landslide mapping 
efforts initiated in Year 1. With the bulk of FOPEN image and DEM analysis following in Years 
2 and 3. 

• North fork of the Gualala River along California’s north coast, for slope erosion analysis in a 
timber harvest site. This data was acquired for work planned to follow in Year 2. 

• The abandoned Iron Mountain Mine west of the town of Redding needs abatement of 
contaminated runoff. Initial data was acquired in Year 1 for work planned to follow in Year 2. 

• An unforeseen opportunity developed to acquire a fifth TOPSAR terrain data set over the San 
Fernando Valley after the January 1994 Northridge earthquake. The mosaiced IFSAR DEM is 
presently being evaluated for use in a FEMA/DMG seismic hazards mapping effort. 

TOPSAR radar data was processed by JPL for the Antioch Hills and Laurel Quad sites. The data 
were processed into three standard products: a radar-derived DEM with 5m horizontal postings; an 
orthorectified SAR image; and, a correlation map. The processed data were recompiled by Calgis and 
provided to DOC for end user experiments and analysis. Calgis developed software to transform the 
data from its spherical system of radar mapping coordinates used by JPL to the cartographic map 
projections and GIS format used by DOC. 


DOC conducted three end-user experiments, two in the Antioch Hills and one for the Laurel Quad 
test site, with JPL and Calgis providing technical support. DOC’s Division of Mines and Geology (DMG) 
Seismic Hazards Evaluation and Zonation Project (SHEZP) conducted a non-vegetated DEM and slope 
hazards experiment in the grass covered rolling slopes of the Antioch Hills and, initiated a vegetated 
DEM and slope hazards assessment in the heavily forested mountainous Laurel Quad, with analysis of 
FOPEN technology planned for Year 2. The Farmland Mapping and Monitoring Program (FMMP) 
conducted a land use and land cover experiment in the urbanized and farmed region of the Antioch Hills 
that is planned to be continued into Year 2 for a change detection experiment. 

4. Antioch Hills Digital Elevation Model Results - Year l 

The IFSAR-derived DEM for Antioch Hills was found to be have a greater accuracy than the U.S. 
Geological Survey (USGS) 7.5 minute quadrangle based DEMs currently available for California. 

Results indicate that a more detailed and more precise analysis of terrain geometry may be possible 
using a radar-derived DEM than with a public domain USGS DEM. Another benefit that radar-based 
DEM production may demonstrate is increased production and reduced turnaround time compared to 
existing methods, this is being studied in Year 2. 


A highly precise DEM was constructed using soft copy photogrammetry and 1:10,000 scale aerial 
photography acquired for a detailed test area, to serve as ground truth and for analysis of the radar- 
derived DEM. SHEZP and Calgis determined a worst case standard deviation of 4.5m vertical which 
makes the radar-derived DEM suitable for creating maps with contour intervals of 15 to 20m. 
Horizontal ground resolution was found to be a factor of 6 times greater or 36 times the data density 
compared to a standard USGS DEM. These results were obtained according to National Map Accuracy 
Standards and procedures based on using a planimetric reference frame from an end-user’s approach. 
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In order to evaluate the TOPSAR data from a sensor performance perspective, JPL performed a 
precise co-registration of a Calgis constructed soft copy-derived photogrammetric reference DEM to the 
TOPSAR generated DEM. After removal of systematic miss registration errors that resulted from 
incomplete calibration of the TOPSAR instrument, JPL found a better level of RMSE accuracy by about 
lm. This indicates that the TOPSAR instrument was performing as expected and after additional 
calibration, may ultimately permit attaining the target goal of 2 to 3m vertical accuracy with a 5m 
horizontal grid. 

5. Antioch hills Land Use and Land Cover Results - Year l 

The IFSAR-derived ortho-rectified image for Antioch Hills was found to have an improved spatial 
accuracy over that of commercially available satellite imagery. A multi-sensor approach of combining 
the IFSAR data with spectral satellite imagery showed potential for improved mapping of agricultural 
land use change and urbanization in California. 


The usefulness of radar imagery as an enhancement to multi-spectral satellite imagery was 
demonstrated in a comparison of four intensity, hue and saturation (IHS) enhanced images developed 
by FMMP and Calgis. Identification of urban, dry vegetation, and marsh land categories was improved 
by substituting an average of Thematic Mapper (TM) and radar principle component (PCI) for intensity 
in an IHS to red-green-blue (RGB) conversion, compared with using TM intensity alone. Discernment 
amongst dry vegetation categories was 25% closer to true ground control with the combined radar/TM 
compared to TM alone. However, the radar/TM combination effected a 15% decline in the ability to 
discern irrigated agricultural categories compared to TM. Enhancement techniques to reduce inherent 
noise and over saturation effects in the radar magnitude image are being pursued as part of Year 2 
efforts. 

Precision point placement tests indicate that the Antioch Hills radar imagery had the highest 
spatial accuracy compared to commercial satellite imagery, and would be suitable for FMMP’s 
registration requirements for mapping farmland and urbanization. The 5m radar pixel makes it possible 
to discern FMMP’s 4 hectare minimum mapping unit which must be measurable in 0.8 hectare 
increments for mapping and reporting agricultural and urban land use changes. 


A five class automated land use and cover classifier (AutoMAP) was developed and demonstrated by 
JPL. Comparison between control land use and the automated classifier results was inconclusive 
concerning which band combinations were the most accurate for automated land use mapping. The 
inherent differences between the land use classes employed by FMMP and the JPL automated classifier 
classes nessicitate further work which is being investigated by FMMP in Year 2. 

6. Year 2 - Overview of Activities 

Year 2 has an expanded scope of work in conjunction with verifying a broader range of applications 
while avoiding duplication of effort. JPL is focusing on development of FOPEN techniques for mapping 
true ground surface using improvements to the TOPSAR instrument and the DC-8 aircraft inertial 
navigation system. In addition, polarimetric data is being acquired and provided using the JPL AIRSAR 
instrument. Calgis is continuing to provide ground truth, photogrammetric and remote sensing 
assistance, along with development of a AIRSAR to GIS convertor. 

Efforts continue with Antioch Hills and Laurel Quad, work on Iron Mountain has been initiated. 
Application expansion and avoiding duplication necessitated discontinuing work on Gualala River. In 
addition to DOC’s internal projects, end-user application has been expanded to include other state and 
local government agencies, and university researchers. 


123 



To broaden the scope of application testing and verification, seven new sites involving ten 
cooperators have been added under the oversight of DOC. Under AKPA’s direction, data from the 
Environmental Research Institute of Michigan (ERIM) X-band IFSAR-E and the Naval Research 
Laboratory (NRL) P-3 aircraft UHF/VHF sensors has been made available for use in Year 2. 


Table 1. Year 2 Expanded Applications 


Site 

Sensor 

End User 

Purpose 

Camp Creek 

JPL-TOPSAR 

UC Davis, Sierra Natural Ecosystems Project 

forest ecosystem study 

Elk Horn Slough 

JPL-TOPSAR 

Moss Landing Marine Labs, Coastal Commission 

estuarine .ecosystem study 

Big Sur 

JPL-TOPSAR 

California Department of Transportation 

coastal highway slope failure 

Camp Roberts 

ERIM-IFSAR-E 

California Army National Guard 

forest fire and range management 

Guadalupe Dunes 

NRL-P3 

DOC Division of Oil and Gas 

mapping buried oil pipelines 

Upper Santa Clara Rive 

ERIM-IFSAR-E 

Los Angeles County Dept, of Public Works 

flood plain assessment 

Lower Santa Clara Rive 

ERIM-IFSAR-E 

Ventura County Flood Control Department 

flood channel erosion control 

Callegaus Creek 

ERIM-IFSAR-E 

Ventura County Flood Control Department 

hydrologic runoff planning 


Systematic errors due to a byte-slip and multi-path interference were noted in the preliminary 
shipments of Year 2 TOPSAR imagery. JPL has developed new calibrations and is reprocessing the first 
data sets that were delivered in 1995. Changes in delivery schedules have lead to modification of 
planned workflow. The primary components of the project have continued despite slipping the original 
plans. 

As a follow-up to Year 1 DEM investigations in the Antioch Hills, SHEZP will be completing its non- 
foliated terrain comparison of the TOPSAR DEM to conventional, large-scale mapping products, which 
include a 1:40,000 scale photogrammetric DEM, a USGS DEM, and a DEM produced by SHEZP using 
the USGS 7.5 minute quadrangle topographic map. A slope algorithm optimal for use on the TOPSAR 
DEM in non-foliate terrain is also under development. A time, cost, and effectiveness study of the 
various available DEMs compared to using radar for DEM construction is being done. Pending progress, 
SHEZP intends to analyze JPL’s preliminary FOPEN DEM and associated data sets for accuracy and 
suitability for use in slope failure assessment. 

Also as a follow-on to its Year 1 land use mapping work in the Antioch area, FMMP is developing an 
improved automated multi-sensor classification workflow, and is improving on computing the accuracy 
of automated classification maps vis-a-vis ground truth. This analysis involves the rasterization of 
ground truth information and a grid-based comparison expanded to the entire study area. Modifications 
to this workflow allow for an improved assessment of the JPL AutoMap classification results. 
Refinements to accuracy assessment procedures will provide a more complete perspective on the values 
of various image enhancement techniques. FMMP will also be working on orthorectification of 1994 and 
1995 imagery for improved ground truth maps, generation of land use conversion statistics between the 
two dates, and a general review of polarimetric data and other products acquired for the Antioch area. 

The Office of Mine Reclamation has initiated work for Iron Mountain, intensive ground truth has 
been accomplished. Calgis acquired close-range photogrammetric images for two areas approximately 
100m x 200m each at the base of the large waste rock fan. Close-range photogrammetry stereo-models 
will be utilized to produce high resolution DEMs as well as for ground truth for RMS height variations - 
comparisons will be made to the radar roughness signature. It is planned to compare hydrologic 
simulation runs using a USGS DEM as a base versus the FOPEN L-band radar DEM. Comparison will 
be based on adequacy of development of drainage patterns and development of rivers and streams. 
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7. Conclusions 


Results of the first year’s efforts have shown that IFSAR can be employed for terrain and land use 
analysis in an environmental setting. A valuable future topic is the development of automated landslide 
hazards and land use change detection algorithms using radar-derived DEMs and images. There exists 
the possibility to train the computer to recognize landslide features using this much more precise radar 
DEM. Such a capability would significantly expedite the process of landslide inventory, quickly covering 
large areas, identifying candidate landslide features and allowing the analyst to concentrate on making 
final judgment calls. Along similar fines, radar-derived automated land use and cover change detection 
warrants additional investigation. Additional work is necessary to evaluate the economics of cost, time, 
and accuracy of existing systems be compared to using radar for DEM construction, and as a part of a 
multi-sensor approach to land use mapping. 
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